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Resumo 

Defeitos topologicos terao sido formados no universo primitivo sendo as suas propriedades 
dependentes dos detalhes da quebra espontanea de simetria que os gerou. Neste trabalho nos 
consideramos consequencias cosmologicas de paredes do dommio, cordas cosmicas e monopo- 
los magneticos. 

Paredes de dommio, formadas quando uma simetria discreta e quebrada, sao o exemplo 
mais simples de um defeito topologico. Neste trabalho investigamos a possibilidade de redes 
de paredes de dommio poderem explicar a actual expansao acelerada do universo. Discuti- 
mos condicoes necessarias para a obtencao de uma rede de paredes de dommio bidimensional 
frustrada e propomos uma classe de modelos que, no limite de grande numero N de campos 
escalares acoplados, aproxima-se do entao chamado modelo 'ideal' (em termos de seu poten- 
cial para produzir frustracao da rede). Utilizando os resultados das maiores e mais precisas 
simulacoes de teoria de campo tridimensionals de redes de paredes de dommio com juncdes, 
encontramos evidencias incisivas para uma aproximacao gradual a uma solucao invariante de 
escala cujos parametros apresentam uma ligeira dependencia em N. Conjecturamos que, apesar 
de ser possfvel a construcao (a mao) de redes estaveis, nenhuma destas redes seriam produtos 
directos de paredes de dommio formadas em transicoes de fase cosmologicas. 

Cordas Cosmicas e Monopolos Magneticos podem ser formados quando simetrias contmuas 
U(l) e SU(2) sao espontaneamente quebradas, respectivamente. Investigamos cordas cosmi- 
cas e monopolos magneticos em modelos tipo-Bekenstein e mostramos que existe uma classe de 
modelos que ainda permitem as solu§6es classicas de vortice de Nielsen-Olesen e de monopolo 
de 't Hooft-Polyakov. Contudo, em geral as solucoes estaticas de cordas e monopolos em mod- 
elos tipo-Bekenstein sao diferentes das classicas com a energia electromagnetica dentro de seus 
nucleos gerando varia§6es espaciais da constante de estrutura fina na vizinhan§a dos defeitos. 
Consideramos modelos com um fun§ao cinetica generica e mostramos que constrangimentos 
provenientes do Princfpio de Equivalencia impoem limites muito fortes as variacoes de a in- 
duzidas pelas redes de cordas cosmicas em escalas cosmologicas. 
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Finalmente, estudamos a evolucao das varia§6es espaciais da constante de estrutura fina 
induzidas por perturbacdes de densidade nao-lineares. Mostramos que os resultados obtidos 
utilizando o modelo de colapso esferico para uma inomogeneidade de comprimento de onda 
infinito sao inconsistentes com os resultados de um estudo local usando gravidade linearizada e 
argumentamos em favor da segunda aproximacao. Tambem criticamos a sugestao de que o valor 
de a de regioes colapsadas poderia ser significativamente diferente do valor de a do universo 
de fundo, com base nos resultados obtidos. 
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Summary 

Topological defects are expected to form in the early universe and their properties depend 
on the particular details of the spontaneous symmetry breaking that has generated them. In this 
work we consider cosmological consequences of domain walls, cosmic strings and magnetic 
monopoles. 

Domain walls, formed when a discrete symmetry is broken, are the simplest example of a 
topological defect. We investigated the domain wall networks as a possible candidate to explain 
the present accelerated expansion of the universe. We discuss various requirements that any 
stable lattice of frustrated walls must obey and propose a class of models which, in the limit of 
large number N of coupled scalar fields, approaches the so-called 'ideal' model (in terms of its 
potential to lead to network frustration). By using the results of the largest and most accurate 
three-dimensional field theory simulations of domain wall networks with junctions, we find 
compelling evidence for a gradual approach to scaling, with the quantitative scaling parameters 
having only a mild dependence on N. We conjecture that, even though one can build (by hand) 
lattices that would be stable, no such lattices will ever come out of realistic domain wall forming 
cosmological phase transitions. 

Cosmic strings and magnetic monopoles can arise with U(l) and 577(2) spontaneous sym- 
metry breaking, respectively. We consider cosmic strings and magnetic monopoles in Bekenstein- 
type models and show that there is a class of models of this type for which the classical Nielsen- 
Olesen vortex and 't Hooft-Polyakov monopoles are still valid solutions. However, in general 
static string and monopole solutions in Bekenstein-type models strongly depart from the stan- 
dard ones with the electromagnetic energy concentrated inside their cores seeding spatial varia- 
tions of the fine structure constant, a. We consider models with a generic gauge kinetic function 
and show that Equivalence Principle constraints impose tight limits on the allowed variations of 
a induced by string networks on cosmological scales. 

Finally, we study the evolution of the spatial variation of the fine-structure constant induced 
by non-linear density perturbations. We show that the results obtained using the spherical infall 
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model for an infinite wavelength inhomogeneity are inconsistent with the results of a local 
linearized gravity study and we argue in favor of the second approach. We also criticize the 
claim that the value of a inside collapsed regions could be significantly different from the 
background one on the basis of these findings. 
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Resume 

Des defauts topologiques ont du etre a l'origine de l'univers primitif, ses proprietes etant 
dependantes des details de la casse spontanee de symetrie qui en sont a l'origine. Dans ce travail, 
nous considererons des consequences cosmologiques de murs du domaine, cordes cosmiques e 
monopoles magnetiques. 

Les murs du domaine, formes quand une symetrie discrete est cassee, sont l'exemple le 
plus simple d'un defaut topologique. Dans ce travail, nous investiguons la possibilite que des 
reseaux de murs de domaine puissent expliquer l'actuelle expansion acceleree de l'univers. 
Nous discutons des conditions necessaries pour l'obtention d'un reseau de murs de domaine 
bidimensionnel frustre e nous proposons une classe de modeles qui, a la limite de grand numero 
N de champs scalaires couples, s'approchent alors du modele appele ideal (vue sa capacite a 
provoquer la frustration du reseau). En utilisant les plus grands resultats et les simulations 
plus precises de theorie de champs tridimensionnels de reseaux de murs de domaine avec des 
jonctions, nous avons trouve des evidences incisives pour une approximation graduelle a une 
solution invariante d'echelle dont les parametres presentent une legere dependance en N. Nous 
avons suppose que, malgre le fait de la possibilite de la construction (manuelle) de reseaux 
stables, aucuns de ces dits reseaux ne serait le produit direct de murs de domaine formes de 
transitions de phase cosmologiques. 

Des Cordes Cosmiques et des Monopoles Magnetiques sont formes quand des symetries 
continues U(l) et 577(2) sont brisees. Nous avons investigue des cordes cosmiques ainsi que 
des monopole magnetiques dans des modeles de type Bekenstein et nous avons montre qu'il 
y a une classe de modeles qui permettent encore les solutions classiques de vortex de Nielsen 
Olsen et de monopole de 't Hooft Polyakov. Cependant, les solutions statiques des cordes et 
monopoles dans les modeles de type Bekenstein sont generalement differentes des classiques 
avec l'energie electro-magnetique a l'interieur de ses noyaux, generant des variations spatiales 
de la constante de structure fine pres des defauts. Nous avons considere des modeles avec une 
fonction cinetique generique e nous avons montre que des desagrements provenant du Principe 
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d' Equivalence imposent de tres petites limites aux variations de a. 

Finalement, nous etudions revolution des variations spatiales de la constante de structure 
fine induites par des perturbations de densite non lineaires. Nous montrons que les resultats 
obtenus en utilisant le modele de collapse spherique, pour une inhomogeneite de longueur 
d'onde, sont peu conformes avec les resultats d'une etude locale qui utilise la gravite linearisee 
et argumentes en faveur de la deuxieme approximation. Nous critiquons aussi la suggestion 
selon laquelle la valeur de a de regions collapsees pourrait etre significativement differente de 
la valeur de a de l'univers de fond, sur la base des resultats obtenus. 
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Units 

Throughout this thesis, the space-time metric g^ u will have a signature H . Greek alphabet 

indices run over space and time while Latin indices over space only. Except if stated otherwise 
we shall assume fundamental units in which 

h = c = ks = Gm 2 Pl = 1, (1) 

then all quantities can be expressed in terms of energy in GeV (lGeV = 10 9 eV). The conver- 
sion factors are 

1 GeV = 1.60 x l(T 3 er# = 1.16 x 10 13 K = 1.78 x 10~ 24 g, (2) 
1 GeV' 1 = 1.97 x l(T 14 cm = 6.58 x l(T 25 s. (3) 

The Planck time and mass are approximately 

t Pl ~ 5.4 x l(T 44 s, (4) 
m Pl ~ 1.2 x 10 19 GeV (5) 

Astrophysical distances will usually be expressed in parsecs, with 

lpc w 3.1 x 10 18 cm, (6) 

or Mpc (1 Mpc = 10 6 pc). 
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Chapter 1 



Introduction 



1.1 Overview 

In the last decades, there were important developments in the interface between cosmology 
and particle physics. In this scenario, topological defects have a wide range of cosmological 
implications. Among these, it has been proposed that domain wall networks might be applied 
to explain the accelerated expansion of the universe at the present time [9J. In the context of 
the varying fundamental constants theories [40J, it has been claimed that the cosmic strings, 
magnetic monopoles and other compact objects might generate space-time variations of the 
fine- structure constant, a, in their vicinity. In this thesis we will study these issues. 

In this chapter we first review the foundations, successes and problems of the standard 
cosmological model in Sec. 1 1 .21 In Sec. 11.31 we introduce the scalar fields and the topological 
defects. We also study the static one-dimensional soli ton solution and some models of scalar 
fields in dark energy scenario. Finally, in Sec. I1.4l we review the varying fundamental constants 
observational results. 



1.2 The Standard Cosmological Model 

The Standard Cosmological Model is based on the assumption that the universe is homogeneous 
and isotropic on large scales. The most general form of a line element which is invariant under 
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32 Introduction 



spatial rotations and translations is 



ds = g^dx^dx" = dt 2 - a 2 {t) 



dr 2 



1 — k r 2 



+ r 2 (d6 2 + sin 2 8 dip 2 ) 



(1.1) 



where g^ u is the metric tensor of the Friedmann-Robertson-Walker (FRW) geometry, a(t) is the 
scale factor and t is the cosmic time. Note that the value of k determines if the universe is flat 
(k = 0), closed (k = +1) or open (k = —1). The line element (|l.ll) describes a homogeneous 
and isotropic space-time and can be parametrized in terms of the conformal time coordinate r 

as 



ds 2 = a 2 (r) dr 2 - 



dr 2 



1 — k r 2 



r 2 (d6 2 + sin 2 6 dcp' 



(1.2) 



In order to describe the components of the energy of the universe, we consider perfect 
barotropic fluids, i.e., perfect fluids with a definite relation between pressure p and energy den- 
sity p. These are described by the tensor 



T^u = {p + p) Uf, u u - pg^, 



(1.3) 



where «„ = dx^/ds and for example a gas of photons has p = p/3. 



General Relativity connects the evolution of the universe to its energy content through the 
Einstein equations 

(1.4) 



R^v — -g^u R — 8nG T^v, 



where is the Ricci tensor, i.e., 



D f) per Q per _i_ per pp p/3 per 



(1.5) 



and R = R^is the Ricci scalar. Using the line element (11.11) one has 



H 2 



ttG k 



P 2' 

3 a 2 



k 

H = -AttG (p + p) + ^ 

or 

p + 3H(p + p) = 0, 



(1.6) 

(1.7) 
(1-8) 



1.2. The Standard Cosmological Model 33 



where H = a/ a is the Hubble parameter and p and p are the total energy density and pressure 
of the universe. The dot represents the derivative with respect to t. Note that while Eq. (1 1.6b is 
obtained from the 00 component of Eq. (11.4b . Eq. (|1.7b results from a linear combination of the 
ij and 00 components. On the other hand, Eq. (|1.8b is obtained by considering the covariant 
conservation of the energy-momentum tensor, V M T^ v = 0. Defining the critical density at a 
given time as 

we can rewrite Eq. (11.6b as 

a 2 H z 

where is the density parameter Vl = p/p c . According to Eq. (11.10b . for p = p c , p < p c and 
p > p c , the universe is flat (fi = 1), open (fi < 1) and closed (fi > 1), respectively. 



1.2.1 The Energy Content of the Universe 

Let us assume that the total energy density of the universe receives contribution from three 
components: matter, radiation and a cosmological constant component A, i.e., p = p m + p r + Pa- 
Since the more recent observational results indicate that the universe is almost flat J9l, we can 
write 

Q = Q m + Q A + Q r = l, (1.11) 

where fij = p i / ' p c is the density parameter for each species. 

Firstly, the contribution of the non-relativistic species (p = 0) is parameterized by p m , that is 
composed by a cold dark matter (CDM) component plus a baryonic one, Q m = Vt dm + Vt b . The 
term cold dark matter means that this component is non-relativistic and does not emit or absorb 
light. The Wilkinson Microwave Anisotropy Probe (WMAP) three year data JH combined with 
the Supernovae la results lfT0l - lfT2ll indicate that 

h 2 Q° dm tt 0.111, h 2 n° b ^ 0.023, (1.12) 



where the subscript denotes the present value of the corresponding quantity, and h parametrizes 
the uncertainty of the value of the Hubble parameter (the Hubble parameter is given by H Q = 
100hkms~ 1 Mpc~ 1 with 0.70 < h < 0.73). Another well known evidence for cold dark matter 
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comes from the rotation curves of spiral galaxies. On the other hand, further indirect evidence 
stems from Big Bang Nucleosynthesis (BBN) [fT3l . 

Secondly, the contribution of the radiation given by p r may be due to photons (0 7 ), neutrinos 
(£l v ) and relic gravitons (Vt g ), i.e., Vt r = fi 7 + Vt v + Vt g . One has 

h 2 « 2.47 x 1(T 5 , h 2 Q° u w 1.68 x 1(T 5 , /i 2 fij < KT 11 , (1.13) 

where the bound on the abundance relic gravitons is based on the analysis of the integrated 
Sachs-Wolfe contribution [fT4| and the neutrino abundance given above is for three species of 
massless neutrinos. 

Finally, another dark component is present in the universe: the dark energy. This component 
has an equation of state 

p 1 

u = ^<-- (1.14) 
P 3 

and dominates the energy of the universe today. It is responsible for the recent acceleration of 
the universe. The main evidence for the existence of dark energy comes from the Supernovae 
la observations and the Cosmic Microwave Background. According to the WMAP three year 
results combined with the Supernovae la ones [fT0ll - [fT2ll . 

h 2 n° A ^ 0.357. (1.15) 



Assuming that h ~ 0.72, the fractional contribution of the various components of the energy 
density of the universe today are 

n2m~o.24, ng~o.o2, n^~o.74, fi°~8.oxHr 5 . (i.i6) 



1.2.2 Distances in an Expanding Universe 

Let us consider an object emitting electromagnetic radiation with wavelength A e . In an expand- 
ing universe, the wavelength received by an observer A is not equal to A e . Since a > 0, the 
observed wavelength will be larger than A e , i.e., 

Ao = t4^ d- 17 ) 

a(t e ) 
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We define the redshift z as 

* = ^f^. (1.18) 

If the luminous signal was emitted when the universe had a fraction of a/a of its present size, 
an observer today should observe a redshift given by 

l + z = ^, (1.19) 

a 

where ao is the value of the scale factor today. We will take ao = 1, except if stated otherwise. 

In a flat universe light traveling freely from t = to dt moves a comoving distance equal to 
dt/a. Therefore the total comoving distance the light could have traveled can be written as 

f* dt' 
Jo a(V) 

Accordingly, no information could have propagated further than 77 since t = is the beginning 
of time. In other words, two regions separated by a distance greater than r] are not causally 
connected, r\ is then named particle horizon. In some particular cases, it is possible to ex- 
press i] analytically in terms of a. For instance, for matter-dominated and radiation dominated 
universes, we have that r\ oc a 1 / 2 and rj oc a, respectively. 

The comoving distance between a distant emitter and a local observer is given by 

f to dt' 

X(a) = / 4l7V- (1-21) 



t(a) a (t') 



Eq. (|1.21l) can be rewritten as 

- 1 da' 



= y. ps^r (L22) 

where we have used H = a/ a in order to change the integration variable. For a matter- 
dominated flat universe, one has H oc a~ 3 / 2 , which yields 



( \ 2 
X{z) = -77- 



(1.23) 



Note that for small z the comoving distance is given by z/H , whereas it asymptotes to 2/H 
for large z 

Let us now consider an object with physical size /. The angle 9 subtended by this object is 
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related to its distance by 

d A = - e , (1.24) 
which is called angular diameter distance and is valid for small 9. As the angle subtended is 

e = -—, (i.25) 

a Xa 

the angular diameter distance for a flat universe is 

d A = T^—. (1.26) 

1 + 2 

The angular diameter distance is then equal to the comoving distance at low redshift but de- 
creases with redshift at high redshift. One direct consequence is that an object at large redshift 
appears larger than it would appear at an intermediate one. 

Another way of determining distances is to consider the flux from an object of known lumi- 
nosity. Let L s be the absolute luminosity, i.e., the power emitted by the a source in its rest frame. 
If one neglects the expansion of the universe and considers that the source is at the comoving 
coordinate r = r s and the detector is at r = 0, the received energy flux is 

F=-r^T- 2 , d-27) 

4 7T <3q r| 

where it is assumed that the light is emitted at t s and observed at t . 
For an expanding universe there are two corrections to Eq. (11.271) : 

:F =- A 2 TTl I V ( L2g ) 

A-na^rj (1 + z) 2 

The first factor (1 + z) is due to the redshift of photons as they travel from the source to the 
detector while the second one computes the fractiorQ between the number of photons that are 
emitted by the source n?, and received by the detector rf per unit of time, i.e., 

Tl e 

-i = l + z. (1.29) 



'This correction is derived by considering that two photons which are emitted with an interval of time St arrives 
separated by an interval of time 5t (ao /a s ). 
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Therefore the luminosity density is 

dl = alr]{l + zf. (1.30) 

1.2.3 The Evolution of the Universe 

The covariant conservation equations lead to the following relations 

Pm = PmO ( — ) , Pr = Pro( — ) , Ph = PkO- (1-31) 



a / v a 



Substituting these relations in the Friedmann equation (11.61) one gets 



H 2 = Hi 



9-n° ^ 

a / 



(1.32) 



where Vt^ = — &/(ag Hq). Of course, the history of the universe then depends upon the relative 
weight of the various physical components of the universe. 

Let us first consider a matter dominated universe (0° = = 0). For = 0, it has a 
decelerating expansion forever with a(t) oc t 2 ^ 3 whereas for > or < the universe 
will collapse in the future or expand forever in a decelerated way. Note that we have neglected 

and have considered Vl m = 1. 

On the other hand, for fl" ^ one finds a different destiny for the universe. In particular, 
= and £l° r = yield the expression 



a vn° 



1 /3 

sinh Q yJnlH (t-to)^ 



2/3 



(1.33) 



that interpolates between a matter dominated universe expanding in a decelerated way a{t) oc 
£ 2 / 3 and an accelerating expanding one, that is dominated by a cosmological constant term. 

Going back in time, it is possible to describe the different epochs until the big explosion - the 
Big Bang. The epoch at which the energy density in matter equals to that in radiation is called 
matter-radiation equality and has special significance for the growth of large-scale structure 
and for the CMB anisotropics. The reason is that perturbations in the dark matter component 
grow at different rates in the matter and radiation eras. Specifically, for z > z eq the universe 
is radiation-dominated with the scale factor a{t) oc t 1 / 2 , while for z < z eq is matter-dominated 
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until dark energy starts dominating with a(t) oc t 2//3 . 

1.2.4 Cornerstones of the Big Bang Model 

1.2.4.1 The Cosmic Microwave Background 

The Cosmic Microwave Background (CMB), a picture of the universe when it was 3 x 10 5 years 
oldj^l is one of the cornerstones of the standard cosmological model. 

Penzias and Wilson lfl~5l measured the CMB spectrum at a frequency of 4.08 GHz and es- 
timated a temperature 3.5 K. A crucial fact about the history of the universe described by the 
CMB measurements is that the collisions with electrons before last scattering implied that the 
photons were in equilibrium leading to a blackbody spectrum. Indeed, the blackbody nature of 
CMB spectrum has been studied and confirmed for a wide range of frequencies ranging from 
0.6 GHz [US El up to 300 GHz flU. 

1.2.4.2 The Primordial Nucleosynthesis 

Another pillar of the standard cosmological model is the Big Bang Nucleosynthesis. The nu- 
clear reactions took place when the universe was between 0.01 s and 100 s old. As a result, a 
substantial amount of 4 He, D, 3 He and 7 Li were produced. The predicted abundances of these 
four light elements agree with the observations, as long as the baryon to photon ratio is between 
4 x 10~ 10 and 7 x 10~ 10 (corresponding to 0.015 <Vt b h 2 < 0.026) JED- The recent WMAP 
results indicate that fij h 2 = 0.0023 ± 0.0008 [9]. Note that since the present universe is almost 
flat (11.161) . the results from the primordial nucleosynthesis imply that the energy density of the 
universe is mostly constituted by a non-baryonic component. 

1.2.4.3 The Hubble Diagram and Supernovae Observations 

The Hubble diagram is the most direct evidence that the universe is expanding. The spectral 
shifts of 41 galaxies were known in 1923. Among them, 36 presented a systematic redshift. In 
1929, Hubble showed empirically that the universe was expanding. He found that the velocity 
of the galaxies increases linearly with distance [20J. Using the same principle of measuring 
the distance and redshift for distant objects, current observations use other objects with known 
intrinsic brightness. The results from observations of the luminosity distance of high redshift 
2 This corresponds to a redshift z ~ 1100. 
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supernovae indicate that the universe is accelerating. In 1998 two research groups (Riess et 
al. and Perlmutter et al. lfT0l -[[T2 | ) published the first sets of evidence that the universe is 
accelerating. Type la Supernovae are thought to have a common absolute magnitude M and 
consequently they are considered good standard candles. 

The simplest explanation for the present accelerated expansion of the universe is a compo- 
nent of the energy density that remains invariant during the cosmic evolution known as Cosmo- 
logical Constant, A. Although this explanation appears well defined, one problem arises when 
we try to understand it in the context of particle physics. We would expect that the value of its 
energy density to be equal to the energy of the vacuum in quantum field theory. In other words, 
we would expect it has the order of the typical scale of early universe phase transitions which 
even at the scale QCD is pa ~ 1CT 3 GeV 4 . Instead, the value is of the order of the critical 
density at the present day, that is, p\ ~ lCT 3 eV 4 . This difference between the expected and 
observed values is called Cosmological Constant Problem. 

The apparent magnitude m of a source in an expanding universe is related to the logarithm 
of T given in Eq. (11.281) as 

m-M = 51og( -^- + 25 ) , (1.34) 

\Mpc / 

where M is the absolute magnitude. Taking the apparent magnitude for two supernovae (one at 
low redshift and another at high redshift) IfTOll one has 

1992P -> m = 16.08 z = 0.026, (1.35) 
1997ap -»• m = 24.32 z = 0.83. (1.36) 

Since at low redshift (z « 1), the luminosity distance can be written as 

d L (z) ~ (1.37) 

-"0 

we apply Eqs. (11.341) and (11.371) to the 1992P supernovae data and find that the absolute mag- 



nitude and the luminosity distance are, respectively, M = —19.09 and d L ~ 1.16H { 



-i 
o • 



Using Eq. (| 1.301) we find that the luminosity distance for a flat FRW universe (11.11) can be 
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written as 

d L = i±£ f d z 1 (1.38) 

Ho Jo V / E,^(l + 5) 3(1+ ^ ) 

Accordingly, assuming that the universe contains only matter (Q,° m = 1), the luminosity distance 
is d L ~ OMH^ 1 . On the other hand, for fi^ = 0.3 and tt A = 0.7 we have that d L ~ l^i^ 1 . 
In other words, an universe dominated by a cosmological constant shows a better agreement 
with the observational data. Note that since the radiation-dominated period is much smaller 
than the total age of the universe, the integral coming from the region z < 1000 does not affect 
the total integral (11.381) . 

Further evidence for the existence of dark energy can be found by computing the age of the 
universe. It can be found by rewriting Eq. (11.321) as 

dZ (1.39) 



h q (i + z)[n° m (i + zf + fiO(i + z y + n° A - + z) 2 } 1 / 2 ' 

The investigation of the Globular clusters in the Milk Way shows that the age of the oldest 
stellar populations is t s = 13.5 ± 2 Gyr j[2T| while for the globular cluster M4 this value is 
t s = 12.7 ±0.7 Gyr (221133. 

However, assuming that Vt^ = 1 and = in Eq. (11.391) we get that to = 8 — 10 Gyr 
which is smaller than t s . A possible explanation is that the universe today is dominated by a 
cosmological constant. 

This is illustrated in Fig. Il.ll ftaken from O). 



1.2.5 Some Problems of the Standard Cosmological Model 

Although the Big Bang theory successfully explains the Hubble expansion law, the CMB and 
the abundances of the light elements, the standard cosmological model does have its problems, 
and we now mention some of them. 

First, the horizon problem stems from the large-scale homogeneity and isotropy of the uni- 
verse. The CMB comes from regions which, according to the standard cosmological model, 
have never been causally in contact. The temperature of the radiation from two different re- 
gions is the same to within at least one part in 10 4 . Consequently, within the standard Big Bang 
model, the early universe was highly homogeneous and isotropic on scales much greater than 
the particle horizon. 
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Second, the question: Why is the present universe so close to flat? summarizes the flatness 
problem. Indeed, the energy density of the universe today is very close to the critical density. 
On the other hand, the critical density is a point of unstable equilibrium, i.e., any deviations 
from Vt = \ grow in time. Accordingly, if today Q ~ 1 then this total density parameter was 
fine-tuned such that |1 — fi| < 1CT 58 at the Planck time. 

Finally, we emphasize that phase transitions in the early universe could have generated 
topological defects. For instance, the spontaneous symmetry breaking whose vacuum man- 
ifold contain non-contractible two-surfaces gives rise magnetic monopoles. This process is 
always present in Grand Unified Theories (GUT), where the symmetry group SU(3) x U(l) is 
generated independent of the original initial group or the intermediate stages of the symmetry 
breaking Il25ll . The existence of heavy monopoles is an inevitable prediction of GUT theories 
and the predicted abundance of this topological defect in the standard cosmological model is 
not consistent with the observations |[T9ll . This problem is named monopole problem. 
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1.3 Scalar Fields and Topological Defects 

The water freezing or melting is an example of a first order phase transition, where there is 
an associated discontinuous change of an order parameter. A second order phase transition 
has a continuous order parameter which is not differentiable. Like phase transitions in con- 
densed matter systems, cosmological phase transitions are associated with spontaneous symme- 
try breaking [|25l . 

Scalar fields can be used to describe the phenomenon of spontaneous symmetry breaking 
ll26l . Despite their simplicity, they have a key role in the understanding of particle physics, 
condensed matter and cosmology. 

Let if be a real scalar field which depends on the space-time coordinates x u , where v = 
0, 1, 2, 3, described by the action 



S = / d}x \J—g 



(1.40) 



where V(<p) is the known A<^ 4 potential given by 

V^) = \( V 2 -V 2 f (1.41) 

with A > and r\ real parameters. 

The high temperature effective potential] for (11.411) can be written as ll25l 

V eS (<p,T)=m 2 {T)<p 2 + j<p\ (1-42) 

where T is the temperature and m is the effective mass of the field ip in the symmetric state 
<p = given by 

m(T) = A ^t 2 - 6t7 2 . (1.43) 
We define the critical temperature T c = y/Qf] such that 

• For T = Tc, the mass vanishes. 

• For T > Tc one has m 2 (T) > and the minimum of V c s is at ip = 0. Since the 
expectation value of cp is null, the symmetry is restored. 



3 We neglected the (,0-independent terms. 
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Figure 1.2: The \(p 4 potential M.41\ for A = rj = 1. 



• For T < Tc one has m 2 (T) < 0. Since if = becomes an unstable state, the field 
develops a non-zero expectation value, say cp = [(Tq — T 2 )/6] 1 / 2 . Given that cp grows 
continuously from zero as the temperature decreases below T c , this process is defined as 
a second-order phase transition. 

The potential (11.411) describes the most elementary type of topological defects (see Fig. 
11.21) . The vacuum manifold, i.e., the set of minima of the potential forms a discrete set and 
intersection of these domains give rise to domain walls that are solutions of the equation of 
motion 



Uip + Av?(y? 2 - rf) = 0. 



(1.44) 



1.3.1 Static One-Dimensional Soliton Solution 



Let us now focus on the static case. The equation of motion is now 



ip dV(ip) 



(1.45) 



dx 2 



dip 



whose one-dimensional solutions must obey the boundary conditions 



x^±oo dx 








lim <f{x) — > ±77, 



(1.46) 
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in order that the energy 



E 



dx 



(1.47) 



be finite. 

The solutions are 




y = ±77 1 arm | \l —i](x — a;o) 



(1.48) 



which are named kink(+) and antikink{-) and take ip from — 77 (77) at x = — oo to r] (—77) at 
x = 00. These solutions are centered at x an d have width[j given by 5 ~ (-\/A 77) -1 . Fig. 11.31 
shows the kink solution for A = 77 = 1 and x = 0. Note that these solutions are Lorentz 
invariant and can be boosted up to arbitrary velocities. We also stress that there is an important 
topological conservation law ensuring the stability of the soliton solution. The topological 
current is defined as j 7 ' = e^ u d u (p, with the associated conserved charge 



Q 



dx f 



00 



<p(x 



-00 



(1.49) 



Therefore the kink (antikink) solution gives rise to Q 7^ 0, that is named topological charge. 
The one-dimensional solution given by (11.481) can be embedded in (3, 1) dimensions. It is then 
a wall separating two regions with different domains (vacua) in the space. These topological 
defects are named domain walls and carry surface tension which is identified with the energy of 
classical solutions in one-dimensional space. 
Let us consider that V(ip) can be written as 



dip 



(1.50) 



tioij] 



where W(<p) is another function^ of <p. In this case the energy is 



(1.51) 



where 



E 



B 



\W[p{x -> 00)] - W[ip{x 



-00 



(1.52) 



4 The width of the topological defect determines the region in which the solution deviates from the vacua. 
5 In supersymmetric scenario, the function W((f) is named superpotential. 
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Figure 1.3: The one-dimensional soliton solution for the Xip 4 potential given in Fig. \1.2\ (we have 
assumed r\ = A = 1 and xq = 0). 

is the Bogomoln'yi energy bound [27 J. In this case, the solutions (11.481) can be found by solving 
the first order differential equations 

£ = (..53) 

ax dip 

and are named BPS (Bogomoln'yi, Prasad, Sommerfield) states j30l[3Tl. It was shown in Refs. 
[1281 - 11331 that Eqs. (|1.53l) factorize completely the equation of motion (11.451) . 



1.3.2 Scalar Field Models and Dark Energy 

Let us consider a homogeneous universe described by a scalar field ip minimally coupled to 
gravity. The Lagrangian density is given by 

C=^d IM (pd*<p-V{<p), (1.54) 

where V(ip) is the potential. 

The energy momentum tensor of the field ip is found by varying the action (11.401) with 
respect to the metric g^, 

T^ = --L/^, (1.55) 
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where 5^g = -(1/2) y/=g g^J g^ v . It gives 



1 



^ K d <pd K <p + V(<p) 



Hence, the energy density in a flat Friedmann Robertson Walker universe is 



(1.56) 



(1.57) 



while the pressure is given by 



• 2 



(1.58) 



where a dot indicates the derivative with respect to time. 

The acceleration of the universe is obtained by combining Eqs. (11.61) and (|1.7I) . It gives 



a 
a 



8ttG 



[P 2 ~ V(ip)] , 



(1.59) 



which implies the existence of an accelerated universe is only possible if ip 2 < V((f). 
Now we recall that the equation of state for the field is given by 



2V 



ip 2 + 2V 



while the continuity equation is 



P = Poe 



./3(l+^)f 



(1.60) 



(1.61) 



where p is an integration constant. Note that the equation of state can vary between the values 
[—1,1]. If My = —1, one has the slow-roll limit, (p 2 « V(<p) and p = const, while io v = 1 
implies Lp 2 >> V((p), in which case the energy density evolves as a~ 6 . On the other hand, 



— 1/3 is the border of acceleration and deceleration. Note that acceleration is realized if 



the energy density is p oc a m , with < m < 2. 



1.3.3 Inflation 

A period of very rapid expansion in the early universe may lead to a solution of some of the 
problems of the standard cosmological model. This was named inflationary scenario and was 
first discussed by Guth P4l . In the inflationary regime: 
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• Regions initially within the causal horizon are blown up to sizes greater than the present 
Hubble radius. 

• The initial curvature decreases by a huge factor, yielding a universe locally indistinguish- 
able from a flat universe. 

• All topological defects formed before inflation, for instance the magnetic monopoles, are 
diluted by an enormous factor. 

1.4 Varying Fundamental Constants 

Are there fundamental scalar fields in nature? Several decades of accelerator experiments did 
not find them, even though they are a key ingredient in the standard model of particle physics. 
We recall that the Higgs particle is supposed to give mass to all other particles and make the 
theory gauge-invariant. We emphasize that there are a number of ways to change the standard 
model of particle physics in order to introduce a space-time variation on the fine-structure con- 
stant a 11351 l36*ll . Since this constant measures the strength of the electromagnetic interaction, 
there are many different ways in which measurements of a can be made. To name just a few, 
locally one can use atomic clocks [|37| or the Oklo natural nuclear reactor Il38ll39l . On the other 
hand, on astrophysical and cosmological scales a lot of work has been done on measurements 
using quasar absorption systems II401 - II43TI and the Cosmic Microwave Background [|44l - [|48ll . 
These different measurements probe very different environments, and therefore it is not triv- 
ial to compare and relate them. Simply comparing at face value numbers obtained at different 
redshifts, for example, it is at the very least too naive, and in most cases manifestly incorrect. 
Indeed, detailed comparisons can often only be made in a model-dependent way, meaning that 
one has to specify a cosmological model and/or a specific model for the evolution of a as a 
function of redshift. Simply assuming, for example, that a grows linearly with time (so that its 
time derivative is constant) is not satisfactory, as one can easily show that no sensible particle 
physics model will ever yield such a dependence for any significant redshift range. 

In order to explain the accelerated expansion of the universe, efforts have been made to 
elaborate models in which this acceleration is related to the variation of the fundamental con- 
stants. In particular, one can reconstruct the dark energy equation of state from variations in 
the fine-structure constant for a class of models where the quintessence field is non-minimally 
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coupled to the electromagnetic field. It has been claimed that variations on a would need to 
be measured to within an accuracy of at least 5 x 10~ 7 to obtain a reconstructed equation of 
state with less than a twenty percent deviations from the true equation of state for < z < 3 



1.4.1 Constraints on variation of a 

1.4.1.1 The Astrophysical Bounds to variation of a 

The investigation of space-time variations of the fine-structure constant has dramatically in- 
creased due to the results of Webb et. al J4TJ. The optical spectra of quasars are rich in absorp- 
tion lines arising from gas clouds along our line of sight. Using the many-multiplet method, that 
exploits the information in many wavelength separations of absorption lines with different rel- 
ativistic contributions to their fine structure, considerable gains in statistical significance were 
achieved. From a data set of 128 objects at redshifts between 0.5 and 3, it is shown in Ref. PTI 
that the absorption spectra were consistent with 

^ se = (-0.57 ± 0.10) xlO" 5 , (1.62) 

where a is the present value of the fine-structure constant. However, the study of 23 absorption 
systems from VLT-UVES quasars at 0.4 < z < 2.3 developed by Chand et. al. 11431 found a 
result consistent with no variation of a, i.e., 

— = (-0.6 ±0.6) x 10~ 6 . (1.63) 

a 

This result was obtained by using a simplified version of the many-multiplet method and con- 
cerns remain about calibrations and the noisiness of the data fit. Also, Quast et. al. and Lek- 
shakov et. al. 11541 observed single quasar absorption systems. They found results consistent 
with 

— = (-0.1 ± 1.7) x 10" 6 (1.64) 
a 



and 



— = (2.4 ±3.8) x 10~ 6 (1.65) 
a 



at z = 1.15 and z = 1.839, respectively. 
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1.4.1.2 Constraints on variation of a from CMB 

The value of the fine- structure constant at high redshift can be measured by investigating the 
CMB. Since a is expected to be a monotonous function of time 11551 . at this epoch any variations 
relative to the present-day value are expected to be larger compared to low redshift observations. 
In fact, a variation of a would alter the ionization history of the universe and thus would have 
an impact on the CMB p4ll49l . Since the redshift of recombination would be changed due to 
a shift in the energy levels of hydrogen. Also, the Thomson scattering cross-section would also 
be modified. 

According to Refs. Il44l - ll47l and ll50l . the overall limit on the variation of a is 

— < 10- 2 , (1.66) 

a 

for z ~ 10 3 . This bound it is consistent with the limit found by Big Bang Nucleosynthesis 
(z ~ 10 9 — 10 10 ) ll56l I571l . It is expected that the Planck mission will provide significant 
improvements in these bounds. 



1.4.1.3 The Oklo natural reactor bounds 

Oklo is the name of the place of a uranium mine in Gabon, West Africa, that supplied the 
uranium ore to the French government. In 1972, the abundance of 235 U was found to be some- 
what below the world standard. This indicated that a self-sustained fission reaction took place 
naturally in Oklo about 1.8 Gyr ago, during the period of Proterozoic, well before Fermi has 
invented the artificial reactor in 1942. Actually, natural reactors had been predicted by Paul 
K. Kuroda, 17 years before the Oklo phenomenon was discovered [|58l . Investigations of this 
phenomenon led to the conclusion that the resonant capture cross-section for thermal neutrons 
by 149 Sm at z ~ 0.15 had created a 149 Sm/ 147 Sm ratio at the reactor site that is depleted by the 
capture process 

149 Sm + n^ 150 Sm + 7 (1.67) 

to an observed value of only about 0.02 compared to the value of about 0.9 obtained in normal 
samples of Samarium. The need for this capture resonance to be in place two billion years ago 
at an energy level within about 90MeV of its current value leads to very strong bounds on all 
interaction coupling constants that contribute to the energy level. This result was noticed by 
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Results of recent investigations using realistic models of natural nuclear reactors provide 
strong bounds on the time evolution of a. Gould et. al use the epithermal spectral indices as 
a criteria for selecting such realistic reactor models [60]. Using numerical simulations, they 
calculated the change in the 149 Sm effective neutron capture cross section as a function of a 
possible shift in the energy of the 93.7 meV resonance. They found that a possible time variation 
in a over 2 billion years should be inside the range 

Aa 

-0.24xl(T 7 < <0.11xl(T 7 . (1.68) 

a 

On the other hand, Petrov et. al used recent version of Monte Carlo codes (MCU REA and 
MNCP) for constructing a computer model of the Oklo reactor taking into account all details 
of design and composition [f6TTl . They obtain conservative limits on the time variation of the 
fine-structure constant: 

Aa 

-4xl(T 8 < < 3 x 1(T 8 . (1.69) 

a 

1.4.1.4 Meteoritic limits on the variation of a 

Peebles and Dicke studied the effects of varying a on the /5-decay lifetime Il62l by considering 
the ratio of Rhenium to Osmium in meteorites, 



lt 7 Re Os + z7 e + e - (1.70) 



which is very sensitive to the value of the fine-structure constant. The analysis of the new 
meteoritic data together with laboratory measurements of the decay rates of long-lived beta 
isotopes indicated a time averaged limit of 

— = (8 ±16) x 1(T 7 (1.71) 

a 

for a sample that spans the age of the solar system (z < 0.45) (631 • 



Chapter 2 

Domain Wall Networks and Dark Energy 



2.1 Overview 

Topological defects are expected to form in the early universe ll64l . Their properties depend on 
the particular details of the spontaneous symmetry breaking that generated them. The choice 
of the minimum after the phase transition is uncorrelated on distances greater than the horizon, 
and defects are necessarily formed. Domain walls, formed when a discrete symmetry is bro- 
ken, are the simplest example. Their small scale properties, such as their physical width and 
tension are fixed and completely determined by the field theory model and are independent of 
the cosmology. 

The average energy density of a frozen domain wall network is given by 

a A a 2 i 
p= — (x — = a , (2.1) 

V a A 

where a is the energy per area unit, A is the wall surface area and V is the volume occupied by 
the domain wall network. 

Hence the equation 

^ + 3Hp(l + w)=0 (2.2) 
at 

implies that the equation of state of a frozen domain wall gas[] is w = —2/3. More generally, 
we have to take into account the dependence of the equation of state of domain walls on their 



'For cosmic strings p oc a 2 , which yields w = —1/3. 
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velocity v. One finds It65ll 

w = ^ = -\ + v 2 . (2.3) 
P 3 

In this chapter we will consider a domain wall network as a possible candidate to explain 
the accelerated expansion of the universe at the present time. We have shown in Chapter 1 
that observational evidence from supernova observations and CMB indicates that the universe 
is accelerating today, i.e., 

w < ji+^m, , _ 1/2 (24) 

Therefore, if a domain wall network frustrates, i.e., is frozen in comoving coordinates (v ~ 0) 
it might explain the present accelerated expansion of the universe [65 J. 

The simplest domain wall networks described by a single scalar field were extensively stud- 
ied and the results indicate that they reach a scaling regime ll66l [6*71 . For this reason, in this 
thesis we consider the possibility that the dynamics of more complex domain wall networks 
produced by two or more coupled fields can lead to frustration. We focus on the role of junc- 
tions, and we study several models where different kinds of junctions can appear. Specifically, 
we consider a number of realizations of the two classes of models introduced in ||69l[6*8l . The 
similarities and differences between them will allow us to point out the key mechanisms at play, 
and further characterize the differences between models with stable Y-type junctions, models 
with stable X-type junctions, and models where both types can co-exist. The details of the 
numerical simulations that we shall present are analogous to those in ll66l [6*71 ITUl . We then 
go beyond these particular models, and investigate an ideal class of models from the point of 
view of obtaining frustrated networks. We study the properties of such networks and investigate 
some realizations by simulating two and three-dimensional networks. 

This chapter is organized as follows. In Sec. I2.2l we present the Press Ryden Spergel (PRS) 
algorithm code and introduce an analytic model for describing the domain wall evolution. In 
Sec. 12.31 we study some models of domain wall network with junctions. In particular, we 
consider the BBL and Kubotani model and present the results of the simulations carried out 
in matter dominated era for two-dimensional boxes. In Sec. 12.41 we investigate the wall lattice 
properties considering the geometry, energy and topology of the domain wall network generated 
by random initial conditions in order to understand what the characteristics of an ideal model 
from the point of view of frustration. We also consider the stability of square and hexagonal 
networks constructed by hand. In Sec. 12.51 we introduce a realization of the ideal model and 
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present the snapshots of numerical simulations for some particular cases. Also, we study the 
scaling behavior of the ideal model and compare it with the BBL model. Finally, in Sec. I2.6l we 
summarize our results and discuss our no-frustration conjecture. 



2.2 Domain Wall Evolution 

Consider the Lagrangian density given by 

C=^d lt <p& t <p-V(( P ), (2.5) 

where V(ip) is the potential. Variation of the action with respect to ip in a Friedmann-Robertson- 
Walker universe leads to the equation of motion 

d 2 cp dip 9 dV 
where V is the Laplacian in physical coordinates and H is the Hubble parameter. 



2.2.1 PRS code 

We simulate the domain wall network evolution using the algorithm introduced by Press, Ryden 
and Spergel in Ref. 117011 in the comoving thickness of the domain wall is fixed. This modifica- 
tion to the equation of motion was introduced in order to be able to resolve the domain walls in 
a comoving grid. 

The equation of motion (12.61) in comoving coordinates can be written as 

d 2 ip / d\na\ 1 dip „ 9 dV 

+ 2 — — - - VV = -a -z-, (2-7) 



dr 2 \ d In r / r dr dtp 
where r is the conformal time. This equation can be generalized to 

d 2 p> „ f d\na\ 1 dip ^ £2 dV 



dr 2 ~^ ^ l \ d In r / r dr ^ ^ dip' ( ^ 

with constant c\ and c 2 . The choices C\ — 3 and c 2 = ensure energy-momentum conservation 
and constant comoving thickness, respectively. 
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Figure 2.1: The matter-era evolution of the domain wall network with of \(p 4 potential (Eq. 17.471 ) and 
a parameter choice A = rj = 1. The simulation starts with random initial conditions inside the range 
(±1). From left to right and top to bottom, the horizon is approximately 1/16, 1/8, 1/4 and 1/2 of the box 
size respectively. 

We integrate Eq. (|2.8I) using finite-difference methods on two-dimensional and three-dimens- 
ional boxes. We assume the initial value of (p to be a random variable between ±t] and the initial 
value of dip/dr to be equal to zero everywhere. We normalize the numerical simulations so that 
r\ = 1. In addition, we set the conformal time at the start of the simulation and the comoving 
spacing between the mesh points to be respectively T{ = 1 and Ax = 1. Fig. 12.11 shows the 
snapshots of a two-dimensional 512 2 simulation carried out for the A<^ 4 potential (|1.41l) . Note 
that the characteristic scale of the network is roughly proportional to the horizon. 

2.2.2 An Analytic Model for the Domain Wall Evolution 

Let us consider an analytic one scale model to describe the domain wall evolution Il66ll67l . 

Firstly, we notice that for a planar domain wall the momentum per unit comoving area is 
proportional to a -1 . Since the physical area in a given comoving volume increases proportion- 
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ally to a 2 , one has jv oc a 3 . Therefore 

^L=(l-tf)(-3Hv). (2.9) 

Secondly, if we consider a network of planar domain walls and the interactions between 
the domain walls are neglected the average energy density of such domain wall network is 
proportional to 7a -1 . Thus it is straightforward to write 

^ = -(l + 3v 2 )(Hp). (2.10) 

Finally, we define the characteristic length L of the network as L = a/ p. Using Eq. 12.101 
we obtain 

^ = (l + 3v 2 )(HL). (2.11) 



2.2.2.1 A More Realistic Scenario 

We have introduced in Refs. JH-Q a more realistic description of domain wall evolution that 
takes into account the curvature of the domain walls, energy loss mechanisms and interactions 
with other fields. 

Firstly, since the wall is not planar there is an acceleration term inversely proportional to 
the curvature radius of the walls, which we assume to be equal to L in the one scale model. 
Therefore the equation for the wall velocity (12.91) can be written as 

* = (1-^(^-31/.), (2,2, 

where k w is the curvature parameter. 

Secondly, we take into account that beyond the damping caused by the Hubble expansion, 
there is a further damping term coming from friction due to particle scattering off the domain 
walls. Hence we define the damping length scale id as 

y = 3H + £ f , (2.13) 

where If accounts the frictional effects due to the interactions of the walls. 
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The equation for the wall velocity (12.121) is now 



< lr r 1 _j ) (k»_Vi (214) 



dt V L 



' d 



where £ d includes both the effects of Hubble damping and particle scattering} Note that the 
friction term if is equal to 

£ f = oc a 4 , (2.15) 

; N W T A v 

where T is background temperature and N w is the number of light particles changing their mass 
across the walls, [64J. Therefore if is expected to dominate at early times, while the Hubble 
term dominates at late times. 

Thirdly, since in a realistic scenario the walls interact there are energy losses to radiatiorj] 
Eqs. (12.1 II) thus become 



dL L r. 

— = HL + -v 2 + c w v, (2.16) 
dt i d 



where the wall energy loss parameter c w parametrizes such energy loss. 

Consider c w = and if -> oo. If a oc t a with a > 1/4, Eqs. (12.141) and (12.161) admit a 
linear scaling solution, with L oc t and v = constant. For example, deep into the radiation 
dominated epoch (a = 1/2) one has 

L = -^=k w t, v = -^=, (2.17) 

whereas deep into the matter era (a = 2/3) one has 

L = ^-k w t } v=-j=. (2.18) 

Moreover, since no strong signatures on the cosmic microwave background have been ob- 
served, the dark energy should be approximately homogeneous and isotropic on large scales. 
The CMB temperature fluctuations generated by domain walls have to be smaller than 1CT 5 
down to scales of the order of 1CT 2 H~ l . This implies that L < lOkpc which leads to a very 



2 If there are no particles scattering off the walls, If — > oo, only the damping from the Hubble expansion is 
computed. 

3 We expect the probability of a wall of size L encountering another of the same size within a time dt to be 
proportional iovdt/L ll64l . 
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small curvature parameter k w < 1CT 6 . This clearly shows that the simplest domain wall scenario 
without junctions, with k w ~ 1, cannot be the dark energy, regardless of any other considera- 
tions. Note that c w ^ leads to larger values of L and does not help to frustration, i.e., L oc a 
and v = 0. Furthermore, the addition of the friction term also does not help much due to the 
limited amount of energy with which domain walls can interact conserving energy and momen- 
tum ll66l [67l . Hence the only possible candidates for dark energy are non-standard networks, 
that is those with junctions where walls intersect. 



2.3 Domain Wall Networks with Junctions 

In the previous section we have shown that only domain wall networks with junctions can 
be considered as candidates to describe dark energy. We have studied several models where 
different types of junctions are formed. In this section, we will investigate two cases: the 
Bazeia-Brito-Losano and the Kubotani models Il68ll69ll . 



2.3.1 The BBL Model 

The Bazeia-Brito-Losano (henceforth referred as BBL model) model was introduced as a pos- 
sible explanation of why the universe presents only three spatial dimensions [1731 - 117711 . The 
Lagrangian density in the BBL model is given by [68] 



1 n 

i=i 

where cpi are real scalar fields and n is a integer number. The potential is 

1 n ( 2 \ 2 
V =2^V~~) +1^1,^2,...,^), (2.20) 
i=i ^ r ' 

where we define the coupling function C between the fields as 

C(<pu ¥2, -, ¥n) =C(cpi, if 2 ) + C(<pu <Ps) + ■■■+ C((pi-u <pi) + ...C(yv_i, cp n ), (2.21) 



with 



C{(fj, if k ) = ~ (y) + ¥ 2 k - Zpfal + (2.22) 
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Figure 2.2: Depending on whether the thick wall has a tension smaller/larger than twice that of the lower 
tension ones the formation of a Y-type/X-type junction will be favored on energetic grounds. 

where k = 1, 2, n. Of course, if the constant coupling e vanishes, one has a system 
formed by n decoupled fields. 

2.3.1.1 BBL Model - 1 field 

The BBL model with just one field is given by 




(2.23) 



with minima ip = ±r (Z 2 symmetry). The potential (|2.23l) is similar to the A(/> 4 given by Eq. 
(Eq. [L4TT) with r = y/2/X = y^\/2r] 2 . (See Fig. |2H> . 

2.3.1.2 BBL Model - 2 fields 

Consider the two-field BBL model given by 

Vfa) = \ £ (r - ^ ) 2 + \ K + <p\ - + 9) . (2.24) 

The minima of potential (12.241) are displayed at the vertexes of a square in the plane (cpi, <^ 2 ) 
whose orientation is defined by the parameter e. 

Firstly, let us assume e = 0. In this case, the fields are decoupled and the potential (12.241) is 
reduced to 

i=l ^ ' 

whose minima cpi — cp 2 — ±r. Therefore the vacua are located at the vertexes of a square 
in the plane cpi,(f2 (Za symmetry). They form six possible topological sectors connecting 
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Figure 2.3: Types of junctions (bottom) and configuration of the minima (top), as a function of the 
parameter er 2 in the two-field BBL model / |6#1/ . 



two different minima: four edges (connecting two neighboring minima in field plane) and two 
diagonals (linking two opposite minima). 

The one-dimensional equations of motions are ll68l 

— = r (2.26) 

ax r 

with the soliton solutions (i = 1, 2) ip\ a = r 2 and tp\ a = r 2 tanh(x) 2 for the edges sectors, and 
Pit — f2 b and ip1 a = r 2 tanh(x) 2 for the diagonal ones, where a = 1, 4 and b = 5, 6. 

Secondly, let us now consider e ^ 0. The potential is well defined only for — 2 < er 2 < 1, 
where if —1/2 < er 2 < 1 the minima are 



^ 2 = T— -J, i = l,2, (2.27) 
1 — er" 4 



r 



2 



while for —2 < er 2 < —1/2 the minima are 



^ = TT^75' *k = °- (2 ' 28) 



Except if stated otherwise we assume r = a/3/2 throughout. The stability of the junctions 
is determined by the choice of e because it determines the tension of the walls. In Ref. |[68l the 
authors find that the tension of the one-dimensional soliton solutions in the limit of very small 
e are E a = 2 + 21s and E b = 4 + 6e, for the edges and diagonal walls respectively f78ll . 

Although they have selected the specific case where only stable X-type junctions survive we 
are interested in all possibles cases: stable Y-type and X-type junctions. 

We confirmed their results in the small e limit, where the ratio of the tension of the walls of 
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Figure 2.4: The matter-era evolution of the domain wall network for the BBL model \2.24\ with r = 
a/3/2 and e = 0.2. The simulation starts with random initial conditions inside the square whose vertexes 



are the vacua (±yl5/7, ±\/l5/7). Note that stable Y -type junctions are preferred. From left to right 
and top to bottom, the horizon is approximately 1/16, 1/8, 1/4 and 1/2 of the box size respectively. 



edge and diagonal types, o e and a d respectively, is given by 

Orf _ 2 + 3e 
a e 1 + 21e/8 ' 



(2.29) 



joining the vacua at ipf = (2/3 — e) -1 / 2 [68J. Depending on whether the diagonal walls have 
a tension smaller (e > 0) or larger (e < 0) than twice that of the edge ones, the formation of a 
Y-type or X-type junction is favored on energetic grounds, as it is illustrated in Fig. 12.21 

Although the analytic expression (12.291) is valid only for e — > 0, the energetic argument is 
valid generically. In the general case, 

• X-type junctions: a d > 2a e for —1 < er 2 < 0. 

• Y-type junctions: a d < 2a e for er 2 > and er 2 < —1. 



Fig. 12.31 illustrates the complete region allowed for er 2 and the respective subsets that yield Y- 
type and X-type junctions. Nevertheless, if the angle between the walls is close to ir/2 then X- 
type junctions may still be formed for v2 o e < o d < 2 a e . Hence, there will be an intermediate 
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Figure 2.5: Same as Fig. \2.4\ except that now e = —0.2. Here the simulation starts with initial conditions 
chosen randomly in the square whose vertexes are the vacua (±yl5/13, i-v/15/13). Note that in this 
case only stable X-type junctions are present in the network. This simulation can be interpreted as the 
two-field Kubotani model given in Eq. \2.49\ with the parameters A = 0.3, rj = y/5/3 and £ = 5/6. 



case for which Y-type junctions are favored but X-type junctions will still be present (see Figs, 
and 1277b. 

Now, we stress that er 2 = —1 is not a marginal case between the formation of the X-type 
and Y-type junctions. Actually, the two fields are decoupled as in the case described by the 
potential in Eq. I2.25l (for er 2 = 0). As a result the fields evolve independently. To be clearer we 
can define a new couple of fields ifti and ^ 2 given by 

Vi = ((fx + V2) , (2.30) 
^2 = -j= Oi - <f2) , (2.31) 

and write the potential in terms of these. One finds that these rotated fields explicitly decouple 
for er 2 = — 1. For an interesting example of a marginal case where both types of junctions are 
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Figure 2.6: Same as Fig. \2.4\ except that now e = —0.4. In this case we start with initial conditions 
chosen randomly in the square whose vertexes are the vacua (rb-v/15/7, 0) and (0, ±yl5/7). Here 
again only stable X -type junctions are present in the network. 

allowed, we can consider the model described by a complex scalar field $ 

£ = d^d^ - V($) (2.32) 

with the potential 

V($) = A \<S> N - l| 2 , (2.33) 
where A is a real parameter. This model has vacua at 

$ = e*#, n=0,l,...,JV-l. (2.34) 

We can alternatively write this field as 

$ = |$|e i0 (2.35) 

and it is then obvious that we will have N minima, evenly spaced around the phase <\>. The case 
N = 2 produces standard domain walls and N = 3 produces Y-type junctions, but the case 
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Figure 2.7: Same as Fig. \2.4\ except that now e = —0.8. Here we start with initial conditions chosen 
randomly in the square whose vertexes are the vacua (rt-v/15/4, 0) and (0, ±yl5/4). Note that in this 
case Y -type junctions are favored in the network. 




Figure 2.8: The configuration space distribution for the four timesteps of the simulation in Fig. 
\2.4\ The vertexes of the square in the plane (^1,(^2) are the minima of the potential, given by 
(±yi577),±VW7)- 
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Figure 2.9: The configuration space distribution for the four timesteps of the simulation in Fig. 
2.51 The vertexes of the square in the plane (^1,^2) are tne minima of the potential, given by 
(±715713)^0^713). 



iV = 4 is slightly more subtle. The potential (12.331) has supersymmetric properties, and hence 
the energy of a specific solution depends only on the initial and final vacua. In other words, the 
possible ways of connecting two opposite vacua (directly or through the intermediate vacuum) 
will have the same energy. Hence this case is an example of the scenario where o d = 2 cr e , which 
means that the formation of Y-type junctions is as likely as the formation of X-type ones, and 
the two will always co-exist. This can also be confirmed numerically, and visually the results 
can be approximately described as the superposition of those ones shown in the two above cases 
with e > and e < 0. 

We confirmed numerically all the properties of the two-field BBL model by performing 256 2 
simulations with r = \/3/2. In this particular case: 

• For -| < e < § the vacua are at { ± (| - e)~^ ±(| - e)~^} 

• For-f <e< -\ the vacua are at {0, ±(§ + f )~ 1/2 } and { ± (| + 0}. 

The outcome of two such simulations is illustrated by the four snapshots displayed in Figs. 
12.41 and [231 where we have respectively taken e = ±0.2. Note that in the former case Y-type 
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Figure 2.10: The configuration space distribution for the four timesteps of the simulation in Fig. 



The vertexes of the square in the plane (ipi, <f2) are the minima of the potential, given by (±\/l5/7), 0) 
and (0^^15/7). 



junctions are preferred while in the latter case only X-type junctions survive. We also verified 
that the results confirm that for any chosen value of e > (that is in the range < e < 2/3, 
given our choice of r), the behavior of the system is similar to that in Fig. I2.4t the formation 
of Y-type junctions is always preferred. On the other hand, we confirm our expectation that 
two distinct cases appear for e < 0. For the interval — 2/3 < e < 0, only X-type junctions are 
formed, as can be seen in Figs. 12.51 and [2761 for e = —0.2 and e = —0.4 respectively. However, 
for e < — 2/3 X-type junctions are no longer stable, and instead Y-type junctions are preferred, 
as illustrated in Fig. I2.7l for the case e = —0.8. 

We also consider the behavior of the configuration space formed by the two scalar fields, 
if 1 and if 2- For this purpose we plot in Figs. 12.8142.1 ll field configurations corresponding to 
the snapshots of Figs. 12 .4142 .7 1 respectively. In agreement with our description, in the cases 
e > and e < —2/3 the field values tend to concentrate along the solutions of both the edges 
and the diagonal sectors of the square. In addition, we stress that the orientation of the square 
is different in the two cases, corresponding to the solutions given by Eq. (12.271) and (|2.28l) 
respectively. On the contrary, for the range — 2/3 < e < the diagonal sector of the square 
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Figure 2.11: The configuration space distribution for the four timesteps of the simulation in Fig. 



The vertexes of the square in the plane (tpi, (^2) are the minima of the potential, given by (±^15/4), 0) 
and (0,±y/15/A). 




Figure 2.12: The configuration of the vacua in the three fields BBL model [68], for the branch A2.37\) 
(Kubotani model H69V for £ > 0) and for (12. 38i (Kubotani model [69] for £ < 0). They are displayed 
at the vertexes of a cube and an octahedron, respectively. Note that in the former case, the formation 
of Y-type junctions is contingent on the relation between the tensions of the walls in the diagonal and 
edge sectors whereas in the latter one, any X-type junctions or Y-type junctions involving the axis sector 
correspond to field space configurations where one of the three fields vanishes. 
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gradually becomes depopulated and only the edges behave as attractors. 
2.3.1.3 BBL Model - 3 fields 

We now move on to the case of models with three scalar fields. The discussion turns out to be 
quite similar to the two-field case, though at some points the effect of the increased dimension- 
ality provides for different phenomenology. The three-field BBL model (12.201) is given by 

- 3e ((p\(p\ + y\ip\ + vlifi) . (2.36) 

As in the two-field case, there are two branches for the minima (see Fig. 12.121) : 

• For — | < er 2 < | there are 8 vacua located at the vertexes of a cube in the space 
(ifi, <f 2 , ^3) w i m values 

There are twenty-eight topological sectors and three kinds of walls, which for obvious 
reasons we can refer to as edges, external diagonals and internal diagonals. The number 
of different walls of each type is respectively twelve, twelve and four. 

• For — 1 < er 2 < — % there are 6 minima which are located at the vertexes of an octahedron 
in the space ((fx, <p 2 , ip$) where 

1 + er z JT ^ 

There are fifteen topological sectors and two kinds of walls, which we can refer to as 
edges and axes. The number of different walls of each type is respectively twelve and 
three. 

Again the choice of the parameter e determines what type of junctions will be present. For 
the branch er 2 > —2/5, corresponding to the cubic solution of Eq. (12.371) . X-type junctions sur- 
vive only if the junctions involving walls from the diagonal sectors are energetically disfavored. 
In other words, they do not necessarily form unless > 3cr e (for the internal diagonals) and 
a e d > 2a e (for the external ones) which is realized only in the interval —2/5 < er 2 < 0, so in 
this range we do have stable X-type junctions. Conversely, for < er 2 < 1/5 the Y-type junc- 
tions are preferred, and these may involve either of the diagonal sectors and so we effectively 
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Figure 2.13: Types of junctions (bottom) and configuration of the minima (top), as a function of the 
parameter er 2 in the three-field BBL model [68 /. Notice that in the case er 2 > there are effectively 
two kinds of Y-type junctions, involving either the internal diagonal or the external diagonal sector. 

have two types of Y-j unctions. On the other hand, for the branch er 2 < —2/5, corresponding to 
the octahedral solution of Eq. (|2.38l) . both Y-type and X-type junctions can survive. Note that in 
this octahedral branch X-type junctions, as well as any Y-type junctions which involve the axes 
sector (as opposed to only the edge sector), correspond to field space configurations where one 
of the three fields vanishes. The cases e = and e = —2/5 are again not interesting because 
they represent the evolution of the three decoupled fields. A numerical study was performed 
along the same lines as what was done for the two-field case, and its outcome is summarized in 
Fig. 

2.3.2 The Kubotani Model 

In this section we investigate the Kubotani model ||69l 1791 and relate it with the BBL model 
studied in the previous section. The potential is given by 

(JV \ 2 N 

E^-^ +^5> 2 -C 2 ) 2 > (2-39) 
i=i / i=i 

where the parameters £ and A are defined such that £ + A > and £ + iVA > 0. 
Let us start by discussing the decoupled case where 

N 

^)=eE(^-C 2 ) 2 - (2-40) 

8=1 

We clearly see that the fields cpi evolve separately: no junctions are formed. The minima of the 
potential (12.401) are (p 2 = ( 2 . Note that when £ = 0, A > the model presents O(N) symmetry. 
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Figure 2. 14: The matter-era evolution of the domain wall network of a Kubotani-type model with a 
parameter choice A = 3/20, rf = 10/3 and £ = 1/12. This also mimics a three-field BBL model with 
r = -y/3/2 and e = —0.2. Note that only stable X-type junctions survive. From left to right and top to 
bottom, the horizon is approximately 1/16, 1/8, 1/4 and 1/2 of the box size respectively. 



The N-th case has a total of 2^ minima, corresponding to the vertexes of an N-dimensional 
hypercube with ipf = rf/N. 

In general, the potential (12.391 ) for N fields has the minima 

2 Ar/ 2 + £C 2 ■ , AT , 0>m 



for £ > and 



2 A?] 2 + ££ 2 2 n /-o n\ 



for £ < 0. Note that the former has 2^ minima, whose vertexes form an N-dimensional (hyper- 
cube), while the latter has 2N minima, that can be thought of as lying in the center of the faces 
of the said N-dimensional cube. 
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Figure 2.15: The matter-era evolution of the domain wall network of a Kubotani-type model with a 
parameter choice A = 3/10, rj 2 = 5/3 and £ = —1/6. This also mimics a three-field BBL model with 
r = and e = —0.4. Note that both Y-type and X-type junctions survive. From left to right and top 

to bottom, the horizon is approximately 1/16, 1/8, 1/4 and 1/2 of the box size respectively. 



2.3.2.1 Relating Kubotani and BBL Models 



Although the Kubotani model (12.391) and the BBL model (12.201) have been introduced in differ- 
ent contexts, they can be related by introducing the generalized potential 

V = Aj2d + BY,tf + Cj2tftf + D, (2.43) 

i i i^j 

where A, B, C and D are real parameters. It is straightforward to show that for 



A = 
C = 



(2.44) 
(2.45) 
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we recover the BBL model while for 

A = A + £, B = -2Ary 2 , (2.46) 
C = 2A, D = At/ 4 , (2.47) 

leads to the potential for the Kubotani model (£ = 0). Therefore assuming 

e = -t = r= vWrm' (2 - 48) 

one relates both models, with the following restrictions 

• For A > and rf > 0, only the interval e < of the BBL model can be mapped into the 
Kubotani model. 

• The assumption £ > is related to the range — 2/5 < er 2 < while £ < is related to 

-1/2 < er 2 < -2/5. 

• The two consistency conditions: 2r/A 2 = 1 on the Kubotani side, and er 2 = — 27e 2 /2 — 
2/3 on the BBL side (though the latter affects only the constant term in the potential) 
indicate that even in this restricted e < range the correspondence is not one-to-one. 

2.3.2.2 Kubotani Model - 2 fields 

The two-field case for the Kubotani model is given by 

V = (A + ( V t + <p*) - 2 r/ 2 A ( V l + (pi) + 2 (pi (pi + X r/ 4 , (2.49) 

that is the same that the BBL model with —1/2 < e < 0. Since the value of the tensions of 
the walls always obeys > 2 a e , only stable X-type junctions are formed here. Note that the 
snapshots of the numerical simulation results shown in Figs. (12.512.91) are the same as for the 
two-field Kubotani model with A = 0.3, rj = a/5/3 and £ = 5/6. 

2.3.2.3 Kubotani Model - 3 fields 

Let us now consider £ = and = 3. Again, the key difference from the BBL models is that 
the tension of the internal diagonal walls is always higher than 3 times the value of the edge 
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Figure 2.16: The matter-era evolution of the domain wall network for a three-field BBL model with 
r = and e = 0.1. Note that only stable Y-type junctions survive. From left to right and top to 

bottom, the horizon is approximately 1/16, 1/8, 1/4 and 1/2 of the box size respectively. 

walls. In other words, unlike the BBL model, Y-type junctions will never be stable in the cubic 
branch because > 2a e d > 3a e . This is shown in Fig. I2.14l for A = 3/20, rf = 10/3 and 
£=1/12. 

Despite the non-standard definition, it is also possible phenomenologically to assume £ < 
0. Here Y-type junctions will form. But just like in the BBL case, one can also envisage 
having X-type junctions, which can arise from field space configurations where one of the three 
scalar fields vanishes. In other words, the scenario with only stable X-type junctions would be 
possible only if the four vacua lying on a plane are selected. However, no realistic spontaneous 
symmetric breaking mechanism that selects only a subset of the vacua is expected to occur. 
Therefore we expect that the existence of both kinds of junctions in any realistic simulation. 
This fact can be confirmed in Fig. 12. 151 for the parameter choices A = 3/10, rf = 5/3 and 
£ = — 1/6. We also note that even though the majority of the junctions are of Y-type, the fraction 
of X-type junctions remains approximately constant during the evolution. For completeness we 
present the case for < er 2 < 1/2 in the BBL model, which has no correspondence in the 
Kubotani case. We have chosen r = \/3/2 and e = 0.1 to illustrate this scenario in Fig. 12.161 
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Figure 2.17: Example of two Y-type junctions in a Kubotani-type model which decay into a stable X-type 
junction; see the text for a detailed discussion. 




Figure 2.18: Example of a X-type junction in a Kubotani-type model which decays into two Y-type 
junctions; see the text for a detailed discussion. 



For further confirmation, we construct by hancj^ a simple square network containing four 
distinct vacua. First, Fig. 12. 171 shows the evolution of £ < Kubotani model starting with four 
vacua lying on a plane and forming two Y-type junctions: these are destroyed and a new X-type 
junction is formed. Although this simulation was performed for the three-field Kubotani model, 
it is valid in the general iV-field case. Regardless of the N chosen, if the vacua located in one 
of the faces of the iV-hypercube are displayed initially as in Fig. 12.171 the Y-type junctions will 
decay into stable X-type ones. In contrast, in Fig. I2.18l we have a set of four vacua that are not 
on a plane in field space: in this case any X-type junction is unstable, and decays into a pair of 
Y-type ones. Again, this conclusion is valid in the general case with arbitrary N. 



4 Construction of a network by hand means selection of special initial conditions in each point into the grid. 
Each point is set to be in a vacuum. 
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2.4 Wall Lattice Properties 

In this section we shall look at geometrical properties of polyhedra in order to investigate if 
equilibrium flat domain wall configurations can be the natural result of domain wall network 
evolution in two-dimensions. We also aim to understand the role played by the model topology 
and the wall tensions on the evolution of the lattice. 

In the following discussion, we will assume that the energy associated with the junctions is 
negligible which in practice means that they are free to move. This is a reasonable assumption, 
at least for the purposes of the present discussion. If it were not the case, one would have to take 
into account the contribution of the junctions when calculating the equation of state associated 
with the domain wall network, and consequently w = p/p would necessarily be greater than 
—2/3 even for a fully static configuration. Such networks would hardly be compatible with 
observational bounds. 



2.4.1 Geometrical Considerations 

One polyhedron is a solid figure with four or more faces, all of which are polygons. In other 
words, a polyhedron is simply a three-dimension solid which consists of a collection of poly- 
gons, usually joined at their edges. Its number of vertexes (V), faces (F), edges {£) and genual 
g surfaces are related by 

V-S + F = 2-20, (2.50) 

that is known as Poincare Formula. We consider simulations in a square box with periodic 
boundary conditions so that g = I. However the choice of boundary conditions will not affect 
our results. 

Let us start by assuming that the number of edges of each polygon, x, and the number 
of edges, d, meeting at a vertex are fixed. Let us denote the number of polyhedron faces by 
F = N x . The number of polyhedron vertexes is V = N x x/d since each polygon has N x vertexes 
but each one of them is shared with d—1 other polygons. Also the number of polyhedron edges 
is equal to S = N x x/2 since each polygon has N x edges but each one of them is shared with 



5 A topological invariant property of a surface defined as the largest number of non-intersecting simple closed 
curves that can be drawn on the surface without separating it. 
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Figure 2.19: Illustration of a random domain distribution in a planar network taken from the fourth 
snapshot of Fig. \2.15\ Each domain has x edges and each junction is formed by d walls. The average 
values are (d) = 53/17 and (x) = 106/19. 

another polygon. As a result, in this case Eq. (12.501) becomes 



This equation has the following solutions (x = 6,d = 3), (x = 4, d = 4), (x = 3,d = 6). 
These are the well known hexagonal type lattices with odd Y-type junctions, square lattices with 
even X-type junctions and triangular lattices with even V -type junctions in 2 dimensions. 

However, in general we do not expect that all the polygons have the same number of edges. 
Therefore we move to the more interesting case where any d is fixed and x is kept free. In this 
case it is straightforward to show that 



with solutions (x) — 6 if d — 3, (x) — 4 if d — 4, (x) — 3 if d — 6 and (x) — > 2 if d — > oo, 
where ( ) represent the average value. If d is not fixed, Eq. (12.521) can be written as 




(2.51) 



(2.52) 



ZZxN* d-2 



((x) - 2) ((d) - 2) = 4. 



(2.53) 



Note that increasing the average dimensionality of the junctions leads to a smaller value of (x). 
We notice that the domains in a realistic domain wall network will not in general have 
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Figure 2.20: An illustration of the decay of an unstable X-type junction into stable Y-type junctions. The 
decay is energetically favorable since it leads to a reduction of the total length of the walls (if all have 
the same tension ). 

straight edges and consequently two edge domains are possible. However, these domains will 
be unstable and collapse due to the domain wall curvature independently of the number of 
elements meeting at each junction. 



Figure 2.21: An illustration of the collapse of three (top) and four (down) edge domains with Y-type 
junctions. The collapse is energetically favorable since it leads to a reduction of the total length of the 
walls. 

In a two-dimensional domain wall network generated dynamically by a cosmological phase 
transition, the number of edges of the domains can assume several values, with a particular 
distribution. The dimensionality of the junctions in such a network, that is, the number of walls 
meeting at each junction, depends on the structure of the chosen potential. We have verified 
this by using the numerical results in Figs. 12.141 (for £ > 0) and 12. 151 (for £ < 0). While in 
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Figure 2.22: An illustration of two different six-edged polygons with Y-type junctions of walls with same 
tension. Both configurations have the same energy. 

the former case only X-j unctions are present, in the latter one there are both stable Y-type and 
X-type junctions. As a simple illustration, we can take the fourth snapshot in Fig. 12.151 which 
is schematically shown in Fig. 12.191 The average number of faces x of each domain is related 
to the dimensionality d of the junctions by 

((x) - 2) ((d) - 2) = 4. (2.54) 

In Fig. 12191 these values are (d) = 53/17 ~ 3.11 and (x) = 106/19 ~ 5.58. We have also 
verified that these numbers are approximately constant during the evolution of the network. 

Finally, we emphasize that no equilibrium configurations exist with d > 6. This means 
that if we started with a domain wall network with d > 6, unstable two edge domains would 
necessarily be present. 

2.4.2 Energy Considerations 
2.4.2.1 Wall Tensions 

Let us suppose that walls with different tensions are present in the network. This assumption 
adds a further source of instability. In fact, the walls with higher tension will tend to collapse 
thus increasing the dimensionality of the junctions which, in general, leads to the production of 
further unstable two edge domains. In the previous section we have illustrated this process in 
Fig. 12.21 which shows the collapse of two Y-type junctions into one X-type junction occurring 
when the thick wall has a tension larger than twice that of the lower tension ones. We also 
simulate the BBL Model (see Fig. 12. 171 and Fig. 12.181) for both the cases with the X-type and 
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Figure 2.23: An illustration of the collapse of a four-edged polygon in the case where two of the sur- 
rounding domains are on the same vacuum state. The collapse leads to a reduction of the number of 
edges of contiguous domains. 




Figure 2.24: An illustration of the collapse of a four-edged polygon in the case where all the surrounding 
domains are on different vacuum states and the domain walls all have the same tension. Again, the 
collapse leads to a reduction of the number of edges of contiguous domains. 

Y-type being energetically favorable. 

In addition, if all the walls have equal tensions, then higher order junctions are unstable and 
decay to stable Y-type junctions. Fig. 12.201 illustrated this point, since the total length of the 
walls (proportional to the total energy of the wall) is smaller for the Y-type junction. 

2.4.2.2 Stability of Domains 

Let us assume that d is fixed to be equal to 3. This will necessarily be the case if all possible 
domain walls have the same tension (see Fig. 12.201) . In this case it is possible to show using a 
local stability analysis that regular polygons with x < 6 are unstable. This is illustrated in Fig. 
12.211 Note that in both cases, x = 3 (top) and x = 4 (bottom), the total length of walls decreases. 
Consequently, the polygons will tend to collapse thus minimizing their potential energy. On the 
other hand in Fig. 12.221 for x = 6, the length remains constant and both configurations have 
the same energy. Hence, given that d = 3 implies (x) = 6, the only possible equilibrium 
configuration with only F-type junctions is a hexagonal type lattice. Otherwise unstable two, 
three, four and five edge domains will be formed, as well as polygons with x > 6. 
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2.4.3 Topology Considerations 

Finally, we investigate the collapse of a particular domain on the network. Consider a square 
with only Y-type junctions with an initial configuration having four contiguous domains with 
Ei, E 2 , E 3 and £4 edges. Since the model has only four vacua, we can, for example, assume 
that the first two will join; the resulting domain will have E\ + E 2 — 4 edges, while the other two 
will have E 3 — 2 and _E 4 — 2 edges after the collapse of the four edge domain. The production 
of three hexagons as a result of the collapse of a four edge domain would require E\ + E 2 = 10 
and E$ = E4 = 8. Hence, we do not expect that a domain wall network in two dimensions will 
naturally evolve towards a hexagonal lattice from realistic initial conditions. We illustrate this 
in Fig. 12.231 We see that the collapse of unstable domains with two, three, four and five edges 
will always result in a decrease in the number of edges of some of the contiguous domains. Of 
course, if there are at least five vacua in the model then the evolution of the local structures 
of the lattice can be different. In fact as the number of vacua increases, the probability of the 




Figure 2.25: The evolution of a perturbed square lattice with even X-type junctions in a matter-dominated 
universe. The top left panel is the initial configuration. From left to right and top to bottom panels the 
horizon is approximately 1/256, 1/20, 1/10 and 1/5 of the box size respectively. The lattice stabilizes in 
the right bottom panel configuration. 
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Figure 2.26: The evolution of a perturbed hexagonal lattice with odd Y-type junctions in the matter- 
dominated epoch. The top left panel is the initial configuration. From left to right and top to bottom 
panels the horizon is approximately 1/256, 1/10, 1/5 and 2/5 of the box size respectively. The lattice 
stabilizes in the right bottom panel configuration. 

2.4.4 Stability of Networks Constructed By Hand 

The stability of simple examples of triangular, hexagonal and square domain wall lattices in 2D 
was studied in f79l l80ll by looking at their macroscopic elastic properties. The authors assess 
the stability of various lattice configurations and find that some of them are stable. Here, we 
confirm and generalize their results numerically showing that there are certain lattice configu- 
rations which are stable against large deformations. Fig. 12.251 shows a field theory numerical 
simulation in the matter epoch for a square lattice with even-type X junctions. We see that 
the network evolves to the minimum energy configuration and stabilizes there. The case of 
the hexagonal lattice with odd Y-type junctions in the matter epoch is shown in Fig. 12.261 As 
was mentioned this is the only possible equilibrium configuration with only Y-type junctions. 
Although hexagonal lattices with Y-type junctions allow for locally confined energy conserving 
deformations, Hubble damping may prevent the collapse of such configurations (even if they 
are perturbed) and consequently one should not completely discard them on this basis. 

Hence, although we agree with the claim that in Minkowski space such a lattice is unstable 
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11791 we find that in an expanding universe it may be stable since domain wall velocities are 
damped by the expansion of the universe. However, we do not expect that hexagonal lattices 
will be an attractor for the evolution of domain wall network simulations. Also, we recall 
that we aim to provide an understanding of domain wall networks generated by spontaneous 
symmetry breaking, and not only the stability of such special configurations. 

These results clearly indicate that the crucial question is not the existence of stable lattice 
configurations but whether these can be the natural result of domain wall evolution with realistic 
initial conditions. 



2.5 The Ideal Model for Frustrated Networks 

In the previous sections we have studied the dynamics of domain wall networks with junctions. 
Here we consider an 'ideal' class of models from point of view of frustration [5 J. This model, 
presented in Ref. [6J, was constructed by considering the energetic, geometrical and topological 
considerations of planar domain wall networks discussed in the previous section. 

The above energy, geometrical and topological considerations led us to propose a class of 
models such that satisfy the following conditions J5l[6l: 

• There are a very large number of minima. In this case the collapse of a single domain 
will only very rarely lead to the fusion of two of the surrounding domains. Note that this 
fusion would lead to a further reduction of the total number of edges of the contiguous 
domains which as we saw is undesirable since it would increase the probability that some 
of the contiguous domains would themselves become unstable to collapse. 

• All possible domain walls have equal tensions. If that were not the case we would be 
adding a further source of instability since the walls with higher tension would tend to 
collapse, thereby increasing the dimensionality of the junctions which, in turn, would lead 
to a decrease in the average number of edges of individual domains and to the production 
of further unstable two-edge domains. 
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2.5.1 Realizations of the Ideal Model 

Let us consider the model given by 

\ N+1 

""lE-fW-S)' (2.55) 

3=1 

with 

N 

^(W-Pij) 2 , (2-56) 



r ■ 

3 

i=l 

where pi. are the N + 1 coordinates of the vacua of the potential. There are iV real scalar 
fields, and the energetic cost for a specific transition between any two of them is the same. 
Geometrically speaking, the p^ are located at the vertexes of a (N + l)-dimensional regular 
polygon. The distance between the vacua is fixed to be equal to r and for large iV we approach 
the ideal model (N — > oo is the ideal case, actually). From the field theory point of view, 
the potential (|2.55l) is the sum of (N + 1) yj 6 -type potentials. Each of them has one minimum 
located at the center, and a continuum of minima at a distance r from the center. Note that N 
of these vacua are located exactly at the centers of the other potentials. 



2.5.1.1 N = 2: the Z 3 Model 

Let us consider the simplest case N = 2 whose potential is written as 



with 



^ = tX>J(''J-''o 2 ) 2 (2.57) 



rf = fa-iy + y*, (2.58) 



Vl + 2j + [ V2 + ^r 1 ^ (2 - 60) 

There are three equidistant vacua at the vertexes of an equilateral triangle with coordinates in 
the {(p x , (p 2 ) plane given by (1, 0), (-1/2, y/3/2) and (-1/2, -y/3/2). (See Fig. [2T271) . This 
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particular case is in fact analogous to the model given by Eq. (12.331) for N = 3, i.e., the Z 3 
model. Fig. 12.281 shows possible collapses of a domain with four edges in the model with only 
three vacua. 



2.5.1.2 N = 3: The Tetrahedral Model 

Let us now consider N = 3. The potential is given by 



^ = tX>J(<'J-<'o 2 ) 2 (2.61) 

3=1 



with 



r{ = (^-lY + ifi + ifi (2.62) 
rl = (<Pi + I) + \<P2 ~ ^] + <pl (2-63) 

Vi + l) + (^2 + ^) +fl (2-64) 



2 / \ T 2 



3 — 

„2 ,2 



^ + ^ + (<^3 - V2) 2 . (2.65) 

There are four vacua located at the vertexes of a tetrahedron as illustrated in Fig. 12.291 Their 
coordinates in the field space (y?i, y? 2 , V^) are (1, 0, 0), (—1/2, \/3/2, 0), (—1/2, —V3/2, 0) and 
(0, 0, V2). The presence of a fourth vacuum allows for the configuration illustrated in Fig. 12.231 
However, the probability that such situation occurs is small and the potential is not enough to 
evade the collapses of domains cutting the network. Another aspect to be noted is that although 
the two-field BBL model also has four minima, the two models are completely different. While 
here only stable Y-type junctions are formed, in the BBL case stable X-type junctions may 
be formed depending on the parameter e. Indeed, if Y-type junctions are formed in the BBL 
model, different kinds of walls necessarily will join the junction. On the contrary, here only 
Y-type junctions are formed and all of them have the same tension. 



2.5.1.3 N = 4: The Pentahedral Model 

Let us assume iV = 4. In this case the five vacua are displayed at the vertexes of a pentahe- 
dron with coordinates in the field space (cpi, <^ 2 , y?3, ^4) given by (1, 0, 0, 0), (—1/2, v^3/2, 0, 0), 
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Figure 2.28: Possible collapses of a domain with four edges in the model with only three vacua. 
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Figure 2.29: Configuration of the vacua of the model ofEq. H2.55i with N = 3. 



(-1 /2, -y/3/2, 0, 0), (0, 0, y/2, 0) and (0, 0, \/2/4, v^/4). Here the potential is written as 



3=1 



(2.66) 



with 



(V?l - 1) + <^ 2 + (f 3 + (f 4 



1\ 2 / \/3 , 

¥>1 + g J + I ^2 + — ) + <^ 3 + ¥>4, 



r 4 = y> 



2 + + (^3 - V^) 2 + </4 



<Pl + <P2 + I ^3 - — I + I <^4 



(2.67) 
(2.68) 

(2.69) 

(2.70) 

(2.71) 
(2.72) 



Now, we have the further possibility of the configuration described in Fig. 12.241 Meanwhile the 
chance of Fig. l2.23l occurring grows. 

In Ref. |[80l . the author introduced a very similar model with five evenly spaced minima (all 
the walls having same tension). Despite the same geometrical structure in the field space, the 
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Carter potential is constructed by summing y? 4 -type potentials, while the ideal model introduced 
here is of y? 6 -type. Of course, both models allow for the configuration in Fig. 12.241 Carter claims 
that this model will form stable X-type junctions, but our results clearly show that only Y-type 
junctions survive. 

2.5.2 Two-dimensional simulations 

Firstly, two-dimensional simulations were carried out to be contrasted with the BBL and Kub- 
otani model results shown in the previous section. A large number of realizations of the ideal 
model in Eq. (12.551) for several values of iV up to 20 were run in the matter dominated epoch. 
In all these simulations we have confirmed that only stable Y-type junctions are formed. In 
addition, if one imposes X-type junctions in the network by hand they are destroyed quickly. 
However, the stable honeycomb network is not formed by starting with random initial condi- 
tions. As an illustration of the simulation results, we have chosen N = A and N = 20. 




Figure 2.30: The matter-era evolution of the domain wall network for the ideal model with N = 4. 
Note that this model is effectively analogous to that of [80], and that only stable Y-type junctions survive. 
From left to right and top to bottom, the horizon is approximately 1/16, 1/8, 1/4 and 1/2 of the box size 
respectively. 



We plot in Fig. 12.301 the snapshots for the evolution of the Pentahedral model (N = 4). 
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Figure 2.31: The matter-era evolution of the domain wall network for the ideal model with N = 20. 
Note that only stable Y-type junctions survive. From left to right and top to bottom, the horizon is 
approximately 1/16, 1/8, 1/4 and 1/2 of the box size respectively. 

This simulation shows the matter era evolution of the domain wall network. Since this case is 
effectively analogous to that the Carter's potential in [80J, it confirms our claim that only stable 
Y-type junctions are formed instead of X-type one as was suggested by Carter. Fig. 12.311 shows 
the corresponding results for N = 20. 



2.5.3 Three-dimensional parallel simulations 

The second part of our results were obtained on the UK Computational Cosmology Consor- 
tium's COSMOS supercomputer. They constitute the most accurate three-dimensional field 
theory simulations in the literature of domain wall networks with junctions and were presented 
in Ref. Q. We have used the same modified version of the algorithm of Press, Ryden and 
Spergel [17011721 described in the previous sections and the code was parallelized with OpenMP 
directives. We have also optimized the code for the shared memory architecture of COSMOS. 
We will discuss the results from series of different eras simulations of 128 3 , 256 3 and 512 3 
boxes, for the ideal class of models (|2.55l) performed for all values of N between 2 and 20. In 



6 See http://www.damtp.cam.ac.uk/cosmos for additional information. 
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addition, we ran 10 different runs for each box size and number of fields. 

These parallel simulations are memory-limited. We need need approximately DIM 3 * N * 
4.5 * 8/1024 2 MB, where DIM is the box size and N is the number of fields. For N = 20 
and 512 3 boxes, not less than 90 GB were required^. The complete set of numerical results will 
appear soon in Ref. JH. An output box binary file can also be produced at specified timesteps 
which can then be used to generate animations, an example of which is available at 

http://www. damtp. cam.ac. uk/cosmos/viz/movies/evo2_25620_msmpeg. avi. 

2.5.3.1 Wall Observables 

We measure the domain wall velocities using an algorithm analogous to that described in [17711 
which removes the radiated energy from the walls which otherwise would contaminate the esti- 
mate of the velocities. This is clearly an important advantage over previous velocity estimations 
(see for example ifTOll ^. 

We defined a domain wall as the region where V(4>) > q V max (V max being the maximum of 
the potential) and estimated the velocities as 



V(<t>i)>0.2V mBX VY ^ ' 



where a dot represents a derivative with respect to conformal time and 7 = (1 — v 2 )~ 1 / 2 . We 
verified that, for 0.2 < q < 0.6, our results are almost independent of q, as long as the domain 
wall is sufficiently resolved. The chosen threshold value q will thus define two properties of the 
network: the corresponding thickness of the static domain wall 5 and the volume fraction of the 
box with domain walls /. They can be used to estimate the comoving characteristic scale of 
the network 

L 5 

L c = - ~ 7 . (2.74) 

a f 

In order to check whether the domain wall network is evolving according to a linear scaling 
solution we define the scaling exponent A so that L c oc r 1_A , where r is the conformal time. If 
A = we have linear scaling. On the contrary, the solution frustrates for A = 1. The scaling 



7 A 512 3 simulation with N = 3 (requiring about 14.5 Gb of memory) takes about 14 seconds per step on 16 
processors and only 5 seconds per step on 32 processors (as the memory ratio becomes favorable) without box 
output, and a complete run takes just 85 minutes. For larger runs the scalability is good if one keeps the memory 
smaller than 1 GB per processor. The largest simulation we have performed, a 512 3 box with 20 scalar fields and 
box outputs at every timestep, took about 2 hours and 15 minutes on 128 processors. 
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Figure 2.32: The asymptotic values ofv* and L c /t for the ideal class of models with N ranging from 2 
to 20. The error bars represent the standard deviation in an ensemble of 10 simulations. The simulations 
were performed in the matter era for 128 3 boxes. 



exponent can be calculated at two different values of the conformal time T\ and r 2 by using the 
following relation ll66l 

A(Tl ' T2)= \n(n/r 2 ) ' (2 - 75) 
where i? = r L^ 1 . In general, since the simulations are carried out until r becomes equal to 
one half of the box comoving size, the values of T\ and r 2 are chosen in order to consider the 
second half of the dynamical range of the simulations. Note that the earlier part is not included 
to avoid contamination of the results by our particular choice of initial conditions, while beyond 
one half of the comoving size, the effect of the finite box size will contaminate the results. 



2.5.3.2 Numerical Results 

We performed series of simulations of 128 3 , 256 3 and 512 3 boxes, for 2 < iV < 20 in matter 
era. We carried out 10 different simulations for each box size and number of fields. 

In Figs. 12.321 and 12.331 we plot the asymptotic values of t>* and L c /t for the ideal class of 
models with 2 < N < 20 run in 128 3 and 512 3 boxes in the matter dominated era, respectively. 
We have performed 10 simulations for each value of N, and shown the standard deviation with 
appropriate error bars IP7T1I . The differences between the results obtained for boxes of different 
sizes are visible but quite small — in general the results are consistent with one another, given the 
error bars. The way the differences go is also what one expects. Mindful of these differences, 
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Figure 2.33: The asymptotic values ofv* and L c /t for the ideal class of models with N ranging from 2 
to 20. The error bars represent the standard deviation in an ensemble of 10 simulations. The simulations 
were performed in the matter era for 512 3 boxes. 
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Figure 2.34: L c as a function of the conformal time t for the class of ideal models ranging from 2 to 20. 
The simulations were carried out in the matter dominated era. Note that as N increases the asymptotic 
value of L c becomes smaller. 
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Figure 2.35: as a function of the conformal time r for the class of ideal models ranging from 2 
to 20. The simulations were performed in the matter dominated era. Note that there is no significant 
dependence of the velocities with N. 




Figure 2.36: The scaling exponents, A, for all N's between 2 and 20 in matter era for 128 3 (dashed 
line), 256 3 (dotted line) and 512 3 (solid line) boxes. Note that frustration would correspond to A = 1 
and linear scaling to A = 0. Note also that as the box size simulation increases, A gets closer to zero. 
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we can still say that for many of our purposes the 256 3 or even 128 3 boxes (each of which we 
can usually run in about an hour or less) will produce adequate results. A larger iV yields a 
smaller asymptotic value of L c /r. This aspect is highlighted in Fig. 12.341 where the differences 
between the successive N results for L c /r become increasingly smaller. Consequently, the 
results obtained up to iV = 20 can be interpreted to be already close to the iV — > oo limit 
results. On the other hand, there is no significant dependence of the velocities, t>*, with N, 
which is confirmed in Fig. 12.351 Overall we find 

v m at = 0.36 ± 0.02 (2.76) 

where v mat is the asymptotic value of the velocity in the matter era simulations. We compare 
(12.761) with the value obtained in Ref. fl\\ for 3D 512 3 simulations of standard domain walls 

V st 

v st ,mat = 0.39 ±0.02. (2.77) 

Therefore the velocities are slightly smaller in the case of networks with junctions. This is of 
course what we expect, though given the various numerical uncertainties the difference between 
the two is not very significant. A more important comparison can be made with the value we 
would expect if the network had no energy losses. In this case we would expect v ~ 0.41, 
and the difference means that energy losses are still noteworthy (at about the ten percent level), 
despite the existence of the junctions. This is ultimately the physical reason why the frustration 
mechanism can never be realized. 

This difference between the behavior of the two averaged quantities is crucial — although 
changing the number of fields N will change the network's characteristic length (and therefore 
its density), the fact that the velocities do not change is an indication that in some sense the 
'local' dynamics of each individual wall will effectively be the same. Note that for this to be 
the case it is obviously necessary that all the walls have the same tension, so we do not expect 
this result to hold for non-ideal models. A second necessary condition is that the junctions 
themselves must be dynamically unimportant. 

Fig. 12.361 shows the scaling exponents, A, for 2 < iV < 20 for 128 3 (dashed line), 256 3 
(dotted line) 512 3 (solid line) boxes. Again, the error bars represent the standard deviation in a 
set of 10 simulations carried out in the matter dominated era. We see that A is slightly greater 
than zero which indicates that there are small departures a the scaling solution (keep in mind 
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Figure 2.37: The asymptotic values ofv* and L c /t for the ideal class of models with N ranging from 2 
to 20. The simulations were performed in the radiation era for 256 3 boxes. 

that exact linear scaling obviously corresponds to A = 0). Indeed we find 

\ mat = 0.11 ±0.05, (2.78) 

where \ mat is the asymptotic value of the scaling exponent in the matter era simulations of do- 
main walls with junctions. Comparing (|2.78l) with the value obtained in Ref. irTTI for analogous 
simulations of standard domain walls 

Kt,mat = 0.04 ± 0.02 , (2.79) 

one has that the two values are in fact relatively similar; the larger exponents in the case with 
junctions should be attributable (at least in part) to the longer time needed for the relaxation 
to scaling when there are several coupled fields. Actually, A is slightly greater than zero but 
as we increase the box size simulations, A gets closer and closer to zero. Larger simulations 
have a larger dynamical range which means that the network evolves for a longer period of 
time. Hence, it is not surprising that those networks are able to approach closer a linear scaling 
solution. On the other hand, frustration (which would correspond to A = 1) is clearly ruled out. 

Even though simulations of the matter era (that is, with the scale factor evolving as a oc 
t 2 / 3 ) are the most relevant from the point of view of cosmological scenarios involving domain 
walls, it is also interesting to study their evolution in other cosmological epochs. Results of 
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Figure 2.38: L c as a function of the conformal time r for the class of ideal models ranging from 2 to 20. 
The simulations were performed in radiation dominated era. Note that as N increases the asymptotic 
value of L c becomes smaller. 




Figure 2.39: as a function of the conformal time t for the class of ideal models ranging from 2 to 20. 
The simulations run in the radiation dominated era. Note that there is no significant dependence of the 
velocities with N. 
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Figure 2.40: The scaling exponents, X, for all N's between 2 and 20 in the radiation era, for 128 3 
(dashed line) and 256 3 (solid line) boxes. Note that frustration would correspond to A = 1 and linear 
scaling to A = 0. Note also that as the box size simulation increases, A gets closer to zero. 



analogous simulations of the ideal model for the radiation era (a oc t 1 / 2 ) are shown in Figs. 
Figs. 12.371 - I2T401 The scaling exponents and properties are defined and measured as before, and 
the simulation parameters are also similar, the only difference being that in this case we only 
have carried out 128 3 and 256 3 simulations. 

The results are qualitatively identical to the ones in the matter era, although as expected 
there are some quantitative differences in the scaling parameters. We have computed the scaling 
exponent X rad : 

\ rad = 0.10 ± 0.05 , (2.80) 

which is again comparable to the one obtained for analogous simulations of standard domain 
walls [rTTTl 

Kt,rad = 0.04 ± 0.02 . (2.81) 
Likewise, we do not find any significant dependence the velocities v ra d with N. We now obtain 

v rad = 0.45 ± 0.03 (2.82) 

whereas for standard domain walls in the radiation era 



Vst,rad = 0.48 ± 0.02 , 



(2.83) 
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Figure 2.41: The scaling exponents, X, for all N's between 2 and 20, for the 128 3 , 256 3 boxes (dotted 
and solid lines respectively) in a slow -expansion era, a oc t 1 / 5 . The error bars represent the standard 
deviation in an ensemble of 10 simulations. Note that frustration would correspond to A = 1. 



and the value we would expect if the network had no energy losses is now v ~ 0.57 (see Eq. 
12.171) . The only noteworthy difference is that energy losses are comparatively more important in 
this case — at about the twenty percent level, as opposed to ten percent in the matter era. There 
is also an enhancement of the scaling density that is similar to the one in the matter era. 

Figs. 12.411 — [27431 display results of analogous simulations for a cosmological epoch with 
an expansion rate even slower than the radiation era, namely a oc t 1 / 5 . This case might be 
cosmologically relevant in some toy models, for a transient epoch in the very early universe, 
but in any case it is an important numerical test because our analytic modeling leads us to expect 
that for a (xt K and k < 1/4 a linear scaling solution can only exist if there are energy losses. 

We emphasize that an accurate numerical evolution of this case is technically difficult. The 
slower expansion on this case means little dissipation, which in turn implies that there is a 
large thermal bath in the box. Nevertheless, we do find evidence for an approach to a linear 
scaling solution in this case, although the evidence is not as strong and unambiguous as that in 
the matter and radiation eras. The main reason for this is that the large thermal bath present 
will tend to make the relaxation to scaling be significantly slower. This explains why the 128 3 
boxes have much larger A than the 256 3 ones- the latter simulations have a dynamic range (in 
conformal time) that is twice that of the former. This analysis therefore confirms the existence 
of energy losses in these domain wall networks, as well as their dynamical importance. 
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Figure 2.42: The asymptotic values ofv* and L c /t for slow-expansion era simulations of the ideal class 
of models with N ranging from 2 to 20 (top and bottom lines respectively). These come from simulations 
c/128 3 , and the error bars represent the standard deviation in an ensemble of 10 simulations. 




Figure 2.43: The asymptotic values ofv* and L c /t for slow-expansion era simulations of the ideal class 
of models with N ranging from 2 to 20 (top and bottom lines respectively). These come from simulations 
of 256 3 , and the error bars represent the standard deviation in an ensemble of 10 simulations. 
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Note that the exponents shown in Fig. I2.41l are effective parameters, calculated in the final 
part of each simulation. Dividing this dynamic range into smaller bins would lead to exponents 
becoming progressively closer to A = as the simulations evolve. The scaling parameters 
and L c /r shown in Figs. 12.421 and 12.431 on the other hand, are not averaged throughout each 
simulations but are measured at the end of each of them. They can therefore still be compared 
with the values we obtained in the radiation and matter eras. As one would expect, in this 
case we find that the wall networks have much higher densities and velocities (v s i w = 0.81). 
The dependence of the scaling density on the number of fields is still seen but is now quite 
weaker than in the matter or radiation cases. Moreover, the scaling velocity now also has a mild 
dependence on N, which as we saw is not present in the other cases. This different behavior is, 
to a large extent, due to the presence of the larger thermal bath, whose effects on a particular 
field will depend on N. The reason for this is the way we set up the simulations. Our choice 
of initial conditions is such that there is always a similar amount of energy in the box, which 
implies that the average energy per field decreases as we increase N. On the other hand, it is 
also conceivable that our algorithm for measuring velocities becomes relatively less precise for 
ultra-relativistic velocities. 



2.5.3.3 Scaling of the Non-Ideal Models 



In order to compare the scaling behavior of the ideal model with the non-ideal ones, we have 
performed some 3D simulations of the three-field BBL Model (|2.36l) in the parameter range 
where only X-type junctions can be formed. We recall that this scenario is realized when the 
parameters er 2 are chosen between —2/5 and 1/2 as in Eq. (12.371) (this case is analogous to 
the Kubotani model with £ > 0). Fig. 12.441 shows the results for L c /r as a function of the 
conformal time r for the BBL model with X-type junctions (dashed line) and the ideal class of 
models for even 2 < N < 20 . Note that the model with N = 2 does not perform better than 
the BBL model from the point of view of frustration, this only happen for N > 3. However, it 
is clear that L c /r — > constant even for large N. as N — > oo. So even the ideal model fails to 
produce a frustrated domain wall network. 
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Figure 2.44: The plot of L c /t as a function of the conformal time r for the three-field version of the BBL 
model (dashed line) with only X-type junctions and for the ideal model (solid line) for even N between 2 
and 20. The simulations were carried out in 256 3 boxes. Although the ideal model with largest N has a 
smaller L c /t than the BBL model, the result shows that the frustration scenario is not realized. 

2.6 Conclusions 

We have studied the possibility that a domain wall network could account for the the dark 
energy or a significant part of it thus explaining the recent acceleration of the universe. For 
this purpose, we considered necessary conditions for frustration. We presented a one scale 
velocity dependent analytic model which describes some of the main properties of the domain 
wall network evolution. We derived a strong bound on the curvature of the walls, which shows 
that viable candidate networks must be fine-tuned and non-standard. We also discussed various 
requirements that any stable lattice of frustrated walls would need to obey. 

We have studied several domain wall models in two spatial dimensions, and discussed the 
conditions under which various types of defect junctions can exist. In particular, we focused on 
the BBL [68] and Kubotani [|69l models. We have highlighted the situations where only stable 
Y-type junctions or X-type ones are formed. In addition we showed the case where both types 
co-exist. 

Using geometrical, topological and energy arguments we isolated model features which we 
expected to be crucial from of the point of view of frustration. However, our results strongly 
lead us to a no frustration conjecture. In all our numerical simulations, including those of the 
ideal model, we find a dynamical behavior consistent with a convergence towards a scaling 
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solution Il66l l67l ITTTl and no evidence of any behavior that could eventually lead to frustration. 
Indeed, 

the results of an 3D simulations clearly indicate that no frustrated network can be formed dy- 
namically out of random initial conditions that would mimic a realistic phase transition. Notice 
that we are making no claims about the possibility of designing (by hand) a frozen lattice-type 
configuration. We have presented the most compelling evidence to date that domain wall net- 
works can not be the dark energy. 



Chapter 3 

Varying-a Cosmic Strings 



3.1 Overview 

The space-time variations of the so-called "constants" of Nature have attracted much attention 
motivated by models with extra spatial dimensions lf8TTl . For example, in the context of varying- 
a, models such as the one proposed by Bekenstein |[82l . the variation is sourced by a scalar field 
coupling minimally to the metric and to the electromagnetic term. 

We recall that the interest in this type of models increased due to the results coming from 
quasar absorption systems pT fl - ll43l suggesting a cosmological variation of a at low redshifts 
ll83l . Other constraints at low redshift include atomic clocks 071 and meteorites [j84]. At 
high redshifts there are also upper limits to the allowed variations of a coming from either the 
Cosmic Microwave Background ll45l - [|50ll or Big Bang Nucleossynthesis [|56l[57l . 

In this chapter we will consider a second class of topological defects, the cosmic strings, 
in the context of the varying-a theories. In the model studied here, the electromagnetic energy 
localized along a stable string-like core acts as a source for spatial variations of a in the vicinity 
of the string. Although a change in the value of a can be accompanied by a variation of other 
fundamental constants, we will consider a phenomenological approach and neglect possible 
variations of other fundamental constants. 

We generalize the results presented initially in Ref. [ 89 ] for the case of a generic gauge 
kinetic function. We find numerical solutions for the fields and compare them with the standard 
solution with no variation of a and compute the spatial variations of the fine-structure con- 
stant in their vicinity (taking into account the limits imposed by tests to the Weak Equivalence 

101 
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Principle ||90ll ). 

This chapter is organized as follows. In Sec. I3.2I we briefly introduce Bekenstein-type 
models and obtain the equations describing a static string solution. We describe the numerical 
results obtained for a number of possible choices of the gauge kinetic function in Sec. I3.3I 
discussing the possible cosmological implications of string networks of this type in the light of 
equivalence principle constraints in Sec. I3.4I Finally, in Sec. I3.5l we summarize our results. 



3.2 Varying-a in Abelian Field Theories 

In this section we describe Bekenstein-type models with a charged complex scalar field <p un- 
dergoing spontaneous symmetry breaking. 

Let us consider a complex scalar field with a gauged U(l) symmetry. The Lagrangian 
density is given by 

C = D^Dy - F*" - V(<j>), (3.1) 
where V(4>) is the potential and the co variant derivatives are defined as 

D^^d^-ieA^, (3.2) 

where e is the electric charge and is the vector gauge field. The electromagnetic field strength 
is defined as 

Fpa, = d^A u - d v A„. (3.3) 



In fact, the Lagrangian in Eq. (13.11) is invariant under the local transformations 

<j>(x) -> <j>'{x) = e lA{x) <j){x), (3.4) 
A^x) -> A'^x) = A^x) + ^ M A(x), (3.5) 

where A(x) is a scalar function of space-time. This set of transformations obeys the U(l) 
symmetry group. 
Varying the action 



S = <Tx 



(3.6) 
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with respect to the fields <p* and A^, one obtains the equations of motion 

dV 

D^Dy+ — = 0, (3.7) 
d,F^=f = -ieWD't- {D»<f>)*<f>}, (3.8) 

where j" is a conserved current. 

In this chapter we will assume the Mexican hat potential V(4>) given by 

y(0) = ^(0*0-^ 2 ) 2 , (3.9) 

where 77 represents the symmetry breaking scale. 

3.2.1 Bekenstein's theory 

Let us now assume that the electric charge e is a function of space and time coordinates, 

e = e e(x> 1 ), (3.10) 

where e(x M ) is a real scalar field and eo is an arbitrary constant charge. 
The Lagrangian density (13.11) is then generalized to include the field e 



C = D^*Dy - - A F, U - V(<f>) - ^ d^e, (3.11) 



where uj is a real parameter and the the gauge invariant electromagnetic field tensor (13.31) is now 

Fp>=\ [d^eAJ-d^eAJ]. (3.12) 

Introducing an auxiliary gauge field given by a M = eA^, the co variant derivatives (13.21) can be 
written as 

Dfi (ft = {df, - ie Q a^) <\>. (3.13) 
In the same way, defining the tensor 

fiiv = eF fll/ = d^a v - d v a^ (3.14) 
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the Lagrangian density can be rewritten as 



£ = D^D^ - -jLfr r - V{4>) 




(3.15) 



which simplifies the variational problem enormously. 

3.2.2 Generalizing the Bekenstein model 

Now, we generalize the Bekenstein model by assuming that the coupling to the Maxwell term 
is not limited to be an exponential function of the real scalar field. 
Let us consider another real scalar field if defined as 



where we have also introduced the gauge kinetic function Bp{ip) = e 2 , and allowed it to be a 
generic function of the field if. In the Bekenstein model, Bp[ip) = exp(— 2 ipj \fu). 

3.2.3 Equations of Motion 

By varying the action with respect to the complex conjugate of 0, i.e., 0*, one gets 




Zj lne. 



(3.16) 



Therefore the Lagrangian density (13.151) becomes 



£ ={p ii <f>Y(D*<f>)--B F {<p)fr,r-V{<l>) + ~d»ipd»ip, 



(3.17) 



dV 
d0* ' 



(3.18) 



Variation with respect to leads to: 



d v [B F (ip) r\ = j 



(3.19) 



with the current j M defined as 



f = ie [0 (Z^0r-0* (Z^0)]. 



(3.20) 
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Note that for B F equal to an arbitrary constant, Eq. (13 . 1 9b reduces to Eq. (13.81) with an effective 
charge e = e / \JB~f~ ■ 



Finally, variation with respect to ip gives 



1 dB F {<p) , 2 



4 



(3.21) 



3.2.4 The Ansatz 



Searching for soliton solutions in this theory, we assume the ansatz proposed by Nielsen and 
Olesen [|9ll , with <p and a M written as 



cf) = x (r)e l 
a e = a(r) , 



(3.22) 
(3.23) 



where x( r ) an d a(r) are real functions of r and all other components of are set to zero. 



Substituting the ansatz given in (I3.22H3.231 into Eqs. (13.1 81 - 13 .2TI) one gets the equations of 
motion 



r dr \ dr 



n \ 2 rf X A 2 

r- e ° a ) ~ — + 2 X 



x = o, 



d ( 1 d , \ rn \ 9 

— B F -— (ra) + 2e - - e a ) X = 0, 
dr \ r dr J V r / 

_ld_ f&A + l dB F (cp) (\d_ y _ Q 
r dr \ dr J 2 dt^ \r dr 



(3.24) 
(3.25) 
(3.26) 



where we have taken into account that 



f = 2f r9 f 1 



V0 



"1 d 

r dr 



1 2 



ra 



(3.27) 



3.2.5 The Energy Density 



We can also investigate the dependence of the energy density on the radial coordinate r. For 
static strings the stress-energy tensor takes a diagonal form with the energy density of the vortex 
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being p = T ° = g 00 T 00 = T 00 , with 

where v = ar, while the spatial components of the stress-energy tensor are given by Tj = 
diag(p r ,p e ,p z ) with p z = —p. Therefore, the energy density of the vortex, p, is everywhere 
positive, while the longitudinal pressure p z is negative. In fact, this is also one of the defining 
features of standard cosmic strings. 



3.3 Numerical Implementation of the Equations of Motion 

We look for solutions to the coupled non-linear equations for a static straight string. We first 
reduce equations (l3.24H3.26l) to a set of first order differential equations for numerical imple- 
mentation. 

Let us introduce three new variables 



dx 
dr 
dv 

dr 

dip 

dr 



br . 
s . 



(3.29) 
(3.30) 
(3.31) 



Then by substituting (13.29H3.3H ) into equations (13.24H3.26l) . one gets 



di 

dr 

db 
dr 
ds 
dr 



r 
1 

~B~f 



n — eo v 



r 



r] 2 X A 2 

1 Y 

2 2* 



X 



dB 



dip 



f , Q / n - e v i 2 
s b - 2e x 



s 1 dBp ,„ 
- + --^— 6 2 . 

r 2 dip 



(3.32) 

(3.33) 
(3.34) 



One has a set of six ordinary first order differential equations, which requires, at least, six 
boundary conditions. 
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3.3.1 Boundary Conditions 

Far away from the core the field \ must be equal to r\ (we take rj = 1). Consequently, substitut- 
ing this condition in Eq. (13.241) one immediately sees that v(r) = n/e far away from the core 
(r — > +oo). In addition, there are four boundary conditions at the string core. As the phase 9 is 
undefined for r = 0, we impose that \ must vanish at that point. Also, we set v = at the center 
of the string. Otherwise, the magnetic energy density would diverge at r = 0. Normalizing the 
electric charge such that e = e at the origin, we impose that the gauge kinetic function is equal 
to unity at r = 0. Finally, by using the Gauss law to solve equation (13.261) assuming that there 
are no sources of a variation other than the string, one has that s = at the string core. 

In summary, the appropriate boundary conditions are 

limx(r) = 0, lim x(r) = 1, (3.35) 

r— >0 r— +oo 

limf(r) = 0, lim v(r) = — , (3.36) 

1 — >0 r^oo eg 

lim B F (r) = 1, lim s(r) = . (3.37) 

Therefore we have a two point boundary value problem with four conditions at the origin and 
two conditions far from the core. 



3.3.2 Numerical Technique 

In order to solve this problem numerically we used the relaxation method which replaces the 
set of six ordinary differential equations by finite-difference equations on a mesh of points 
covering the range of the integration. This numerical method is very efficient if a good initial 
guess is supplied. In our case, the solutions of the standard Nielsen-Olesen vortex were used to 
generate a good initial guess. We checked that our code reproduces the results for the standard 
Nielsen-Olesen vortex if B F = 1. 

The results that we will describe here are obtained by assuming in the numerical implemen- 
tation that A = 2 and n = 1 for definiteness and use units in which if = 1. We also emphasize 
that although the Principle Equivalence tests constrain dBp/dLf to be very small, we will often 
use larger values in order to better investigate the effects associated with a variability. 
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Figure 3.1: The numerical solution of the scalar field if? = lne as a function of distance, r,from the core 
of string, for the original Bekenstein model. If uj < Q then e — > oo when r — > oo. On the other hand, if 
uo > e — > w/j^n r — > oo. 772e dashed line represents the constant-a theory, which corresponds 
to the limit uj — » oo. 

3.3.3 Exponential coupling 

The general prescription detailed above can be particularized to specific choices of gauge kinetic 
function. First, let us consider the particular case 

B F {<p) = e~% , (3.38) 

where uj is a coupling constant. Actually, this case corresponds to the original Bekenstein model 
described by Eq. (13.151) . 

It is easy to show that in the limit uj — » oo one recovers the Nielsen-Olesen vortex with 
constant e. Although the gauge kinetic function in Eq. (|3.38l) is only well defined for cu > the 
model described by the Lagrangian density in Eq. (13.151) allows for both negative and positive 
values of to. However, note that if uj < the energy density is no longer positive definite. 

In Fig. 13.11 we plot the numerical solution of the scalar field if? = lne as a function of 
distance, r, from the core of string, in the context of the original Bekenstein model. Note that 
if to < then e diverges asymptotically away from the string core. On the other hand if u? > 
then e goes to zero when r — > oo. In the large uj limit the curves for positive and negative 
uj are nearly symmetric approaching the dashed line representing the constant-a model when 
uj — > oo. 
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3.3.4 Polynomial coupling 

Now let us move to a more general case where the gauge kinetic function takes a polynomial 
form 

JV 

B F (<p) = 1.0 + ^2 fa<p\ (3.39) 
1=1 

where fa are dimensionless coupling constants and iV is an integer. We have solved the equa- 
tions of motion assuming the gauge kinetic function given in Eq. (13.391) , and have verified some 
important properties. 

First, we note that if fa — the classical Nielsen-Olesen vortex solution with constant a is 
still a valid solution. This means that there is a class of gauge kinetic functions for which the 
classical static solution is maintained despite the modifications to the model. It can be easily 
verified by substituting the gauge kinetic function (13.391 ) in Eqs. (13 .24113 .261) and taking fa = 0, 
that the standard Nielsen-Olesen vortex solution is recovered. 

Substituting the gauge kinetic function in Eq. (13.261) one gets 

M i rd l) = 7 (|> - D-W^*) + j (EWA**") . (3-40) 

which shows that the transformation fa — > —fa for odd i modifies the sign of the solution of 
ip(r) without changing x or 6 since Si? is kept invariant. We will see that both for fa > and 
Pi < the behavior of ip = In e is similar to that of the original Bekenstein model described by 
the Lagrangian density in Eq. (|3.15l) . with u > 0, in particular in the limit of small |/5i |/large 
uj. In fact a polynomial expansion of the exponential gauge kinetic function of the original 
Bekenstein model has fa = —2/y/u. This relation between the models arises from the fact that 
B F given in Eq. (13.381) can be expanded in a series of powers of <p according to 

(-2Y 

fa = • (3.41) 

On the other hand, as mentioned before, if fa = one recovers the standard result for the 
Nielsen-Olesen vortex with constant-a, irrespective of the chosen values of fa for i > 1. 

We have studied the behavior of the solutions of Eqs. (13 .29143 341) for various values of N 
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Figure 3.2: The numerical solution of the scalar field ip = lne as a function of distance, r, from the 
core of string, for a polynomial gauge kinetic function. Models 0, 1 and 2 are defined by (3\ = (P2 
arbitrary), Pi = —3, = (linear coupling) and Pi = —5, P2 = 10 respectively. 



but for simplicity we shall only consider iV < 2. In particular, we consider 



with two free parameters. In Fig. I3.2l we plot the numerical solution of the scalar field ij) = In e 
as a function of distance, r, from the core of string, for a polynomial gauge kinetic function. 
Models 1 and 2 are defined by [3\ = —3, Pi = (linear coupling) and (5\ = —5, fa = 10 
respectively. Model (dashed line) represents any model with [3\ = and has a = constant. 
Note that the replacement (5\ — > —fii does not modify the solution for ip. In Fig. I3.3l we plot the 
numerical solution of the fields x(r) and b(r) as a function of distance, r, from the core of string, 
for models 0, 1 and 2. We see that the change in b(r) with respect to the standard constant-a 
result is much more dramatic than the change in x(r). In order to verify the modification to 
x(r) in more detail we define the function 



and plot in Fig. 13 .41 the results for the different models. We clearly see that even a small value 
of Pi leads to a modification of the vortex solution with respect to the standard Nielsen-Olesen 
solution. Finally, we have also studied the behavior of the energy density in this model. In 



B F = 1.0 + p l( p + p 2 ip 2 , 



(3.42) 




(3.43) 



3.3. Numerical Implementation of the Equations of Motion 111 



3.5 - 

3 - 

2.5 - 

S 2 - 
-Q 



1.5 - 




r 



Figure 3.3: The numerical solution of the fields x( r ) an d b(r) as a function of distance, r,from the core 
of string, for models 0, 1 and 2. Note that the change in b{r) with respect to the standard constant-a 
result is much more dramatic than the change in x( r )- 




Figure 3.4: Plot of Ci( r ) = l°g(%i/xo) fo r the different polynomial gauge kinetic functions. Note that 
although this was not very visible in Fig. \3.3\ even a small value of f3\ leads to a change of the vortex 
solution. 
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Figure 3.5: The energy density as a function of the distance, r from the string core for models 0, 1 and 
2. The dashed line represents the constant-a model These results clearly show an increase of the energy 
density due to the contribution of the extra field (p. 

fact, since the fine structure constant varies we have a new contribution due to the field ip, to 
the total energy of the topological defect. As has been previously discussed in ref. Il89l the 
contributions to the energy density of the string can be divided into two components. One 
component is localized around the string core (the local string component), the other is related 
to the contribution of the kinetic term associated with the spatial variations of the fine structure 
constant and is not localized in the core of the string. The energy profile of this last contribution 
is analogous to that of a global string, whose energy per unit length diverges asymptotically far 
from the core. In Fig. I3.5l we plot the string energy density as a function of the distance, r from 
the string core for models 0, 1 and 2. The dashed line represents the constant-a model. We 
clearly see an increase of the energy density due to the contribution of the extra field ip. 



3.4 Cosmological Consequences of Varying-a Strings 

In the previous section we studied static solutions of varying-a cosmic strings. Now, let us 
investigate whether such cosmic string networks can induce measurable space-time variations 
of a in a cosmological setting. In order to answer this question we recall that cp satisfies the 
Poisson equation given by 

V 2 <P = \/3f 2 - (3.44) 
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In Eq. (I3.44I) we have assumed for simplicity that the gauge kinetic function is a linear function 
of (p, i.e., 

B F (<p) = l + P<p. (3.45) 
The variation of the fine structure constant away from the string core is then 

2vrr^ = /3/(r)//(r max ) , (3.46) 
ar 

where 

fi(r) = 2tt f p(r')r'dr\ (3.47) 



.7 

and 

I(r) = 9 r \ f f( r Ydr' , (3.48) 

*l~ l <\ r max) Jo 

is a function of r smaller than unity. Here r max represents a cut-off scale which is in a cosmo- 
logical context of the order of the string correlation length. Far away from the string core J(r) 
is a slowly varying function of r which is always smaller than unity. An approximate solution 
for the behavior of the field ip may be obtained by taking J(r) ~ const 

V? ~ In — , (3.49) 



2tt 

where r is an integration constant. Since e = B f (lp)~ 1 / 2 we have 



■In | 


/ r 




V r o 




-1/2 



e ~ 1 - ^ r -) in ) . (3.50) 



47r \ro 

We see that the variation of the fine structure constant away from the string core is proportional 
to the gravitational potential induced by the strings. 



3.4.1 The upper limit to the variation of a 

In order to compute an upper limit on the variation of the fine- structure constant in the vicinity 
of a cosmic string, we first recall that the value of Gfi is constrained to be small (< 1CT 6 — 10~ 7 ) 
in order to avoid conflict with CMB and LSS results [1771192119311941 (or even smaller depending 
on the decaying channels available to the cosmic string network ll95l ). We stress that these 
constraints are for standard local strings described by the Nambu-Goto action. Even though the 
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cosmic strings studied here are non-standard, having a local and a global component, we expect 
the limits on G/j,(r max ) to be similar as those on Gfi for standard local strings. 

On the other hand, the factor (3 appearing in Eq. (I3.50I) is constrained by equivalence prin- 
ciple tests to be \(3\ < 1(T 3 G 1/2 J90l|9i. 

Hence, taking into account that we cannot observe scales larger than the horizon (~ 10 4 Mpc) 
and are unlikely to probe variations of a at a distance much smaller than 1 pc from a cosmic 
string, Eq. (13.501) implies that a conservative overall limit on observable variations of a seeded 
by cosmic strings is 

— < ID"' 2 (3.51) 
a 

which is too small to have any significant cosmological impact. We thus see that, even allowing 
for a large contribution coming from the logarithmic factor in Eq. (13.501) . the spatial variations 
of a induced by such strings are too small to be detectable. 



3.5 Conclusions 

We have investigated cosmic strings in the context of varying-a models considering a generic 
gauge kinetic function. We showed that there is a class of models of this type for which the 
classical Nielsen-Olesen vortex is still a valid solution. We also found that the spatial variations 
of a (sourced by the electromagnetic energy concentrated along the string core) are roughly 
proportional to the gravitational potential induced by the strings which is constrained to be 
small. However, we showed that Equivalence Principle constraints impose tight limits on the 
allowed variations of a on cosmological scales induced by cosmic string networks of this type. 



Chapter 4 

Varying-a Magnetic Monopoles 



4.1 Overview 

In the previous chapter we have studied a model where spatial variations of the fine-structure 
constant were induced by cosmic strings. In this chapter we study another class of topological 
defects, local magnetic monopoles, which can arise in the context of non-Abelian field theories 
with SU (2) spontaneous symmetry breaking Il97l - ll99l . In the context of varying-a models, the 
electromagnetic energy in their core may also induce spatial variations of a 0. 

We recall that the monopole problem of the standard cosmological model can be solved by 
considering, for instance, the inflationary universe ll34l . Inflation introduces a period of expo- 
nential expansion of the universe, and all of the present observed universe arises from a tiny 
region which was initially smaller than the causal horizon. As a result, not only monopoles, 
but all defects formed before the inflationary period, are diluted by an enormous factor. Af- 
ter inflation the universe thermalizes at some temperature and some monopole- antimonopole 
pairs may be produced by thermal fluctuations, but their density is supressed by a very large 
value HlOOl I101II . However, in some models the predicted monopole abundance produced at 
thermalization is acceptably small, but still detectable H102L A number of other theories were 
elaborated in order to solve the monopole problem. Among them, there is the possibility that 
monopoles and anti-monopoles are linked by strings 11031 which can also lead to the annihila- 
tion of magnetic monopoles. 

In this chapter we start by solving the standard problem for a fixed a and recover the clas- 
sical solution presented by 't Hooft and Polyakov ll97l l98l (see for instance Refs. Il25l l26l l99l 
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I104I0 . Then we proceed to investigate other solutions where a varies in space, comparing them 
with the standard one. Finally, we find the upper limits to the variation of a in their vicinity. 

This chapter is organized as follows. In Sec. I4.2I we introduce Bekenstein-type models 
in Yang-Mills theories and obtain the equations of motion for a static magnetic monopole. 
We write the energy density of the monopole and show that the standard electromagnetism is 
recovered outside the core if a is a constant. In. Sec. I4.3I we present the numerical technique 
applied to solve the equations of motion and give results for several choices of the gauge kinetic 
function. In Sec. I4.4l we study the limits imposed by the Equivalence Principle on the allowed 
variations of a in the vicinity of the monopole, as a function of the symmetry breaking scale. 
Finally, in Sec. I4.5L we summarize and discuss the results. 



4.2 Varying-a in Non-Abelian Field Theories 

In this section we describe Bekenstein-type models in the context of non-Abelian Yang-Mills 
theories H105L We shall consider the non-Abelian gauge field theory described by the La- 
grangian density 

C = D fl ® a D fM <& a - -F^F afU/ + 1/($ a ), (4.1) 

where $ a with a = 1, 2, 3 is a triplet whose elements are scalar fields and V($ a ) is the potential. 
The covariant derivatives are 

D^ a = d^ a - ee abc A b ^ c (4.2) 
and the electromagnetic strength tensor has the form 

F; v = d,Al - d v Al + ee abc A b ^ (4.3) 

where A ^ is an isovector, e abc is the Levi-Civita tensor and e is the electric charge. 

Analogously to the Mexican hat potential V{4>) taken for the string case in Eq. (13.91) . we 
assume the potential 

V($ a ) = - ($ a $ a -r? 2 ) 2 , (4.4) 
where A > is the coupling of the scalar self-interaction, and r\ is the vacuum expectation value 
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of the Higgs field. 



4.2.1 Bekenstein-type models 

Now we proceed as in Chapter [3] and assume that the electric charge is a function of the space- 
time coordinates, i.e., e = e e , where e is a real scalar field. In this case, the Lagrangian density 
can be written as 

C = ±(D ll (D" $ a ) - V m - f% r v + \d,pd^, (4.5) 

where a® = eA a ^ is an auxiliary gauge field and ip is a real scalar field. 
The covariant derivatives are now 

D^ a = d^ a + e Q e abc a b ^ c (4.6) 

and the non-Abelian gauge field strength is given by 

ft, = £ F % = d,< ~ d v al + e e abc afc, (4.7) 

where F" is defined in Eq. (14.31) . We recall that B f (lp) = e~ 2 is the gauge kinetic function 
which acts as the effective dielectric permittivity and can phenomenologically be taken as an 
arbitrary function of ip. 

Note that the Lagrangian density in Eq. (|4.5I) is invariant under SU(2) gauge transformations 
of the form 

e abc $ 6 A c , (4.8) 
e a6c a 6 A c + e ^A a , (4.9) 

where A a is a generic isovector. This symmetry is broken down to U(l) because there is a non 
vanishing expectation value of the Higgs field. Thus the two components of the vector field 
develop a mass Mw = V eo, while the mass of the Higgs field is Mh = rj y/2X. 



5<t> a 



8< 
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It is convenient to define the dimensionless ratio 



M, 



H 



^2X 



Mi 



w 



e 



(4.10) 



and to rescale the radial coordinate by M w , so that distance is expressed in units of M^, 1 . 



4.2.2 Equations of Motion 



Varying the action with respect to the adjoint Higgs field $ a1 " one gets 

Dn <3> a = -2 A $ a (<E> b <E> 6 - rj 2 ) . (4.11) 
Variation with respect to leads to 

d u [B F ( V )r^]=r (4.i2) 

with the current j afl defined as 

j"" = e e abc <5> b (D»<5> c ). (4.13) 

We stress that if B F is a constant the standard theory is recovered with an effective electric 
charge e = e / \fB~F~ ■ 



Finally, variation with respect to (p gives 



1 dB F (<p) , 



d ll & t i P =--r—^-f 2 (4.14) 

4 dip 



with/ 2 = f^r^. 
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4.2.3 The ansatz 

We are now interested in static, spherically symmetric, magnetic monopole solutions. Hence 
we take the "hedgehog" ansatz 





[r) 


= X(r) 


< 


J) 


= o, 


< 


[r) 





X 



Xk 

e r 2 



[W(r) - 1], 



(4.15) 
(4.16) 
(4.17) 



where x a are the Cartesian coordinates and r 2 = x k Xk- X(r) and W(r) are dimensionless 
radial functions which minimize the self-energy, i.e., the mass of the monopole 



E 



4:711] f°° J rr 2 fdX\ 2 „ t2 ( 2 r 2 



1-X 



2\2 



+ B t 



dW 
dr 



+ 



[I - w 

2^" 



2\2 



+-i-Y 

2 V dr J 



(4.18) 



where the coordinate r and the functions X and cp have been rescaled as 

r 

r — > , X — > vX, ip — > vcp. 

e v 

So that r, X and are now dimensionless. 

It will prove useful to compute the energy density which is given by 



(4.19) 



P 



HI 

e I 



Bi 



1 fdW 

2 



i n-w 2 \ 2 ' 

2 



1 /dXV /W*Y C 2 ,, ^2 1 



(4.20) 



Substituting the ansatz given in (I4.15I44~T71) into the equations of motion (14.1 1144~T4l) one obtains 



1 d 

r 2 dr 
d 
dr 

1 d 

r 2 dr 



By 



dr 
dW 
dr 

2 dip 
dr 



2 
Bp 

dB F 

dip 



<\i ^ 2W 2 
'\ - X 2 ) H — 



T 

2\ , v 2 



(1 -W z )+X 



X = 0, 
W = 0, 



1 /l 
2 



(4.21) 
(4.22) 

(4.23) 
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4.2.4 The electromagnetic tensor outside the monopole 

Defining 

U» = w\ + eabc $a ( D ^ b )( D ^ c ) ( 4 - 24 ) 

and choosing the gauge where $ a = 5 a3 |$|, i.e., one gets 

f lt » = d ll al-d v al. (4.25) 

By writing = a M , one identifies (I4.25I) with the usual electromagnetic tensor which for 
constant cp satisfies the ordinary Maxwell equations everywhere except in the region where 
X ~ 0. 



4.3 Numerical Implementation of the Equations 

In order to solve the equations of motion numerically, we reduce them to a set of first order 
equations. For that purpose we define the variables 



dX 

dr 
dW 

dr 
dip 

dr 



V. 



u, 



(4.26) 
(4.27) 
(4.28) 



Equations (I4.21H4.23T) then can be written as 



dV 

dr 
dU 

dr 

ds 
dr 



r 2 V 



V. 



2W 2 X 



1 

~B~f 



^sU + ^(W 2 -l)W + X 2 W 

dtp r z 



2 s dB F 

r dip 



U 



1 - w 



2\ 2 



(4.29) 
(4.30) 

(4.31) 



which require at least six boundary conditions. 

At far distances from the core (r — > oo), the Higgs field X(r) falls off to its vacuum value, 
X — 1. Using this boundary condition in Eq. (14.211) . one gets that W(r) must vanish far from 
the core. On the other hand, since the symmetry at the core is not broken, then the Higgs field 
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vanishes at r = and from the regularity of the energy-momentum tensor, W = 1 at r = 0. 
The other two boundary conditions come from the normalization of the electric charge at the 
origin. At the core e = e , which means that the gauge kinetic function is equal to unity. 
Finally, by using the Gauss law to solve the Eq. (|4.23l) without sources of a variation other than 
the monopole, one gets that s must vanish at r = 0. Therefore, we have at all six boundary 
conditions, four of them at origin and the other two far away from the core: 

0, IimW(r) = l, (4.32) 

i — >o 

1, lims(r)=0, (4.33) 

1, lim W(r) = 0. (4.34) 

r— >oo 

As the boundary conditions are at different points of the domain of the functions to be found, 
we use the relaxation numerical method replacing the set of differential equations by finite- 
difference equations on a grid of points that covers the whole range of the integration, as in the 
previous chapter. 



limX(r) 
lim Bp{r) 
lim X(r) 

r— >oc 



4.3.1 't Hooft-Polyakov Standard Solution 

In this section we take Bp = 1. This recovers the standard 't Hooft-Polyakov monopole solution 
ll97l with no variations of the fine- structure constant, i.e., for which ip(r) = constant. 

It proves to be useful to write the equations of motion which are 

dX 

dr 
dW 



- V, (4.35) 

dr 



^ U, (4.36) 

dr 

dV 2V C 2 , 2W 2 X 

-T = -— + ^-X(X 2 -l) + ^^, (4.37) 

dr r 2 r A 

= \(W 2 - 1)W + X 2 W. (4.38) 

dr r 2 



Let us first consider £ = 0, i.e., a massless Higgs field. In this case the equations of motion 
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are 



dX 

dr 
dW 

dr 
dV 

dr 
dU 

dr 



V. 



U. 



IV 2W 2 X 
+ 



;{W 2 - l)W + X 2 W. 



(4.39) 
(4.40) 
(4.41) 
(4.42) 



Substituting Eqs. (14 .3 91 - 14 .421) in Eq. (14.181) . the energy of the monopole can be written as 



E 



Anrj 
eo <- Jo 



dr 



1 / dX 

2 \ ~dr 



1 - W< 



dr J 



dP\, (4.43) 



where we have introduced a new scalar field P(r) as 



P(r) =X(1 -W 2 ). 



{AAA) 



We notice that the equations of motion (l4.39H4.42l) can be completely factorized into the set of 
two first order equations 



dX 

dr 
dW 

dr 



1-W 2 



-wx. 



(4.45) 
(4.46) 



These are precisely equivalent to the equations of motion (I4.39H4.421) . This is the Bogomoln'yi 
method E71 - ir32l . Using Eq. (14.431) to calculate the energy one gets 



E 



Airi] 
e 



P(r -»• oo) - P(r = 0) 



(4.47) 



which from (14.35H4.38~1) gives E = Anri/eo, that is the Bogomoln'yi energy. Note that the 
solutions of equations (l4.45H4.46l) are 



X(r) 
W(r) 



sinh (r) ' 
1 1 



tanh (r) r ' 



(4.48) 
(4.49) 
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Figure 4.1: The numerical solution of the field X(r) as a function of distance, r,from the core of the 't 
Hooft-Polyakov monopole for several values of Q. For £ = the Higgs field is massless and for £ — > oo 
Higgs field is frozen at its vacuum value except at the origin. 




Figure 4.2: The numerical solution of the field W(r) as a function of distance, r,from the core of the 't 
Hooft-Polyakov monopole for three values of £. 
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Figure 4.3: The numerical solution of the energy of the 't Hooft-Polyakov monopole, E(C)> as a function 
of the ratio of the Higgs mass to the vector field one, (. Note that the monopole energy remains finite 
when C — > oo. 

which are the Bogomol'nyi, Prasad, Sommerfield (BPS) solutions Il27lll06ll . In Figs. I4.1l and 
Owe plot X(r) and W(r) given by Eqs. (I4.48I) and ( I4.49I ). 

Let us now consider the opposite limit for (, i.e., £ — > oo which for fixed eo implies A — > oo, 
i.e., the Higgs potential is much larger than its kinetic term. As a result, the Higgs field is forced 
to be frozen at its vacuum value everywhere except at the origin. The equations of motion (14.351 - 
I4.38I) are reduced to 



and X = 0. 

In Figs. I4.1l and l4~2l we plot X{r) and W(r) for £ — > oo. As expected for £ = the Higgs 
field is massless while it is frozen at its vacuum value (except the origin) for £ — > oo. 

The previous values for ( simplify very much the equations of motion. However, we have 
solved them numerically with a very good precision for ( in the interval 10~ 4 < £ < 10 3 . We 
also show (see Fig. I4.3I) that the energy of the monopole increases with ( asymptoting to a finite 
constant for £ — » oo . 

In the next subsections we will find magnetic monopole solutions with varying fine-structure 
constant. We will set the parameter £ = 1. 




(4.50) 
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Figure 4.4: The numerical solution of the scalar field ip = lne as a function of distance, r,from the core 
of the monopole, in the Bekenstein model. Note that ifu < 0, e — > oo when r — > oo. Ora o?/ier /lancf, 
ifu>0 then e — > w/ien r — > oo. 77ie dashed line represents the constant-a theory, which corresponds 
to the limit uj — > oo. 



4.3.2 Bekenstein Model 



We already know from Chapter [3] that for 

# F (p) = e "^, (4.51) 
with cj > 0, one recovers the original Bekenstein model 11821 . Now defining 




(4.52) 



and substituting in (14.511) one gets that Eq. (l4.5l) can be written as 

C = l - (Dp $ a ) (D" d>°) - ^/;, - y ($ a ) + |<9 M ^ & V, (4.53) 

where now uj can take either positive or positive values. 

First we note that in the limit uj — > oo one recovers the 't Hooft-Polyakov monopole de- 
scribed in the previous subsection. In Fig. 14.41 we plot the numerical solution of the scalar field 
ip(r) = lne for several values of uj. One concludes that if uj < then e diverges asymptotically 
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r 



Figure 4.5: The numerical solution of the energy density p as a function r, the distance from the core of 
the monopole. The dashed-line is for uj — > oo. Note that for uj = 2.0, the energy density is non negative, 
while for uj = —0.1 it becomes negative at some distance from the core. 



away from the core of the monopole, and the energy density 



V l e -2, 

e ' 



dW 
dr 



+ 



1 fl-W 2 \ 2 



dr J \ r J 



(4.54) 



is no longer positive definite. However if uj > then e vanishes asymptotically when r — > oo, 
and the energy density is, in this case, positive definite (See Fig. 14.51) . We also note from Fig. 
14.41 that in the large uj limit the curves for positive and negative uj are approximately symmetric 
approaching the dashed line which represents the constant-a model when uj — ► oo. 



Finally, in Fig. I4.6l we plot the numerical solution of the scalar field X(r) and the gauge field 
W(r) as a function of distance, r, from the core of the monopole. The dashed-line represents 
the constant-a solution and the solid line represents the Bekenstein one for uj = 2.0. We verified 
that even in the uj — » limit, the change in X{r) with respect to 't Hooft-Polyakov solution is 
still negligible. 
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Figure 4.6: The numerical solution of the Higgs field X{x) and the gauge field W(x) as functions of 
distance r, from the core of the monopole. Note that the change in W(r) with respect to the 't Hooft- 
Polyakov solution(dashed-line) is more relevant than the change in X(r). 

4.3.3 Polynomial Gauge Kinetic Function 

Let us now consider a class of gauge kinetic functions: 

N 

B F {<p) =1.0 + (4.55) 

where iV is an integer. As ip is dimensionless then $ is dimensionless as well. We have shown 
in Chapter [3] that these models are equivalent to the Bekenstein one (14.511) since 



This relationship between the coupling constants /3, and u has interesting consequences for the 
model given by Eq. (14.551) . 

First we verified that the behavior of ift both for (5\ > and /3i < is similar to that of the 
Bekenstein model with u> > which is recovered in the limit of small |/3i|/large to. We also 
verified that if one takes [3i = one gets the 't Hooft-Polyakov limit, for any $ with i > 1. 
This means that there is a class of gauge kinetic functions for which the classical static solution 
is maintained despite the modifications to the model. This property also is present in the case 
of the cosmic strings as we have mentioned in the previous chapter. 



Pi 




(4.56) 



w l / 2 i\ 
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Figure 4. 7: The numerical solution of the scalar field ip = In e as a function of distance, r, from the 
core of monopole, for a polynomial gauge kinetic function. Models 0, 1 and 2 are defined by f3\ = (P 2 
arbitrary), Pi = — 3, p 2 = (linear coupling) and Pi = — 2, p 2 = 5 respectively. 

Another property can be noticed when one substitutes the gauge kinetic function in the 
equation of motion for ip (14.231) . One gets 



r 2 dr \ dr J 



C 2 



N 



N 



J2(2k - l)(3 2k „i V 2k - 2 + ^(2k)fo k <p 



2k-l 



k=l 



k=l 



(4.57) 



with 



C 2 



dW 



+ 



W 



2\ 2 



(4.58) 



dr J 2 \ r 

Since that C 2 > 0, when $ — > —Pi for odd i one sees by Eq. (14.571) that <p(r) — > —ip{r). 
However, as Bp in Eq. (|4.55l) is kept invariant, X(r) or W(r) do not vary. 

We have found the set of solutions for several values of N. Without loss of generality we 
took N = 2 in Eq. (14.551) to write 



B F = 1.0 + Picp + p 2 ip 2 



(4.59) 



and define the model with Pi = (standard monopole solution), model 1 with Pi = —3, P% = 
(linear coupling) and model 2 with Pi = —2, [3 2 = 5 (quadratic coupling), respectively. Fig. 
14.71 shows the numerical solution of the scalar field ift(r). As a result we have verified that the 
replacement f3\ — > —Pi does not modify the solution for ip as it was expected. 

From the equations of motion (l4.21H4.23l) a change in Bp is directly manifest in a change of 
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Figure 4.8: The numerical solution of the fields X(r) and W(r) as functions of distance, r, from the 
core of monopole, for models 0, 1 and 2. Note that the change in W(r) with respect to the standard 
constant-a result is much more dramatic than the change in H(r). 



0.02 




r 



Figure 4.9: Plot of Li(r) = log(Xj / 'Xq) for the models 0, 1 and 2. One clearly sees that even a small 
value of Pi leads to a different vortex solution from the standard 't Hooft-Polyakov one. 
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Figure 4.10: The energy density as a function of the distance, r from the monopole core for models 0, 
1 and 2. The dashed line represents the constant-a model. Note that there is an increase of the energy 
density due to the contribution of the extra field (p in the case models 1 and 2. 

U, i.e., of W and indirectly results in a change of V, i.e., of X. This can be checked through Fig. 
14.81 where one can see clearly that the change in W(r) with respect to the standard constant-a 
result is much more dramatic than the change in X(r). This is more explicit in Fig. 14.91 where 
we define the function 



with i = 0, 1, 2. We note that even a small value of [3i leads to a modification of the magnetic 
monopole solution with respect to the standard 't Hooft-Polyakov solution. 

We have also studied the behavior of the energy density in these models which are plotted 
in Fig. 14.101 The increase of the energy density for models 1 and 2 when compared to model 
is due to the contribution of the spatial variations of the structure constant. In fact, as in the case 
of cosmic strings, the energy density of the monopole can be divided into two components: one 
that is localized inside the core of the monopole and another related to the contribution of the 
kinetic term associated with the spatial variations of the fine structure constant. 




(4.60) 



4.4 Constraints on Variations of a 



We now look for an overall limit on the spatial variations of the fine- structure constant 
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seeded by magnetic monopoles. We consider a monopole network with the same core size, r c . 
We impose a cut-off by considering that there is a separation of 2 r max to its first neighbor. In 
other words, we set r max /r c ^> 1 so that we can neglect the interactions of each monopole 
with its neighbor and therefore approximate it as an isolated monopole, such as the one in the 
previous sections, with the difference that the space outside the core now extends up to a finite 
distance: the correlation length of the system r max . 

Let us now calculate the spatial variations of the fine-structure constant outside the core. For 
simplicity we assume 

B F {<p) = 1 + p<p, (4.61) 

i.e., the gauge kinetic function is a linear function in ip that satisfies the spherically symmetric 
Poisson equation given by 

V 2 <P=\pf. (4-62) 

Note that (3 is constrained by Equivalence Principle tests to be such that \(3\ < 1CT 3 G 1//2 (see 
ESI |96l). 



Integrating Eq. (l4.61l) from the core up to r max , the cosmological cut-off scale, one gets 



Ik = PI(r)M(r max ) . (4.63) 
In Eq. (l4.63l) we used the mass of the monopole which is given by 

rr 

M(r) = 4tt / p(r')r' 2 dr' (4.64) 

and 



o 



I ^ = J7T^ I f 2 ( r Y 2 dr, (4.65) 

1VJ - Vmax) JO 

which is a slowly varying function of r outside the core always smaller than unity. Thus we can 
take J(r) ~ const and integrate Eq. (14.631) to get 

<p , (4.66) 



47r \r r 

where r is a integration constant which could be identified as the core radius. Using Bp(<p) 
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e 2 one gets 

t ~l-^»-» (I -IV (4.67) 

8tt \r r / 

which means that the variation of the fine structure constant away from the monopole core is 
proportional to the gravitational potential induced by the monopoles. 



4.4.1 GUT Monopoles 

Let us estimate an overall limit for the spatial variation of a outside the core of GUT monopoles. 
In this context the mass of the monopole is of order of 

M(r max )~10 w GeV. (4.68) 

The variation of a from the core up to infinity is 

A« = a(r ^ oo) - q(r ) = g2 _ x < 10 -i3 (4 . 69) 
a cx{ro) 

where we have used a = a /B F ((p), with a = a(r — > oo) and a = a(r ). 

4.4.2 Planck Monopoles 

Proceeding as above for the Planck scale symmetry with 

M(r max ) ~ 10 19 Gev (4.70) 

one obtains an overall limit for the spatial variation of the fine- structure seeded by magnetic 
monopoles 

^ < ID"', (4.71) 
a 

which is still very small even for Planck scale monopoles. 



4.5 Conclusions 



We introduced Bekenstein-type models in the Non-Abelian Yang-Mills theory in order to in- 
vestigate the variations of the fine-structure constant in the vicinity of a static local monopole. 
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We studied various models with constant a and reviewed the standard static 't Hooft-Polyakov 
magnetic monopole solution. Then we considered models with varying-a and confirmed that 
despite the existence of a class of models for which the 't Hooft-Polyakov standard solution 
is still valid. This property had been verified in Chapter [3] in the context of the cosmic strings 
in varying-a theories. However, in general, the solutions of the local defects depart from the 
standard case. We also showed that Equivalence Principle constraints impose tight limits on the 
variations of a induced by magnetic monopoles. This confirms the difficulty to generate sig- 
nificant large-scale spatial variation of the fine structure constant found in the previous section, 
even in the most favorable case where these variations are seeded by magnetic monopoles. 



Chapter 5 



Evolution of a in the Non-Linear Regime 

5.1 Overview 

We have seen in this thesis that a real scalar field coupled to the matter sector could be re- 
sponsible for the variation of a. In fact, given that in the course of the cosmological evolution, 
inhomogeneities grow, become non-linear and decouple from the background evolution, it has 
been argued that the same could happen to the local variations of a. A number of authors have 
studied this problem and have shown using linear theory that the spatial variations of the fine 
structure constant induced by fluctuations in the matter fields are proportional to the gravita- 
tional potential. Although such spatial variations of a are typically very small to be detected 
directly with present day technology, some authors have argued that large variations could occur 
in the vicinity of compact objects with strong gravitational fields [86J-[88J. 

In previous chapters, we have considered two examples of compact objects, cosmic strings 
and magnetic monopoles. We have shown that both kinds of topological defects do not seed 
large variations of a in their vicinity despite the electromagnetic energy concentrated in their 
core. In this chapter, we will study the variations of a induced by the collapse of large scale 
structures. In fact, this has been previously studied in [|86l - [l88l using a simple spherical infall 
model for the evolution of infinite wavelength density perturbations and a particular generaliza- 
tion of the Bekenstein model [85| for the evolution of a. It was found that in the linear regime 
and in the matter era the variation of a would follow the density contrast. Moreover, it was 
also claimed that this approximation was valid in the non-linear regime (meaning turnaround 
and collapse). Here we revisit and question these results. Our work was divided into two parts. 

135 
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Firstly, we confirm the results obtained by Barrow and Mota by using the infinite wavelength 
approximation to study the growth of a uniform spherical matter inhomogeneity with a final 
size r Il86l - lf88l . Secondly, we argue that a local approximation for determining the variation 
of the fine- structure constant outside virialized regions is the correct approach H107L 

An alternative approximation was followed in Ref. H108L This provides a more detailed and 
mathematically-inclined analysis of local variations in physical 'constants', but does confirm 
our results. Our approximation, while much simpler, has the advantage of making explicit the 
reasons why spatial variations have to be small, and why the use of the spherical collapse model 
of a infinite wavelength perturbation is inadequate. Other papers also have later confirmed our 
results ifTO^ -lfTTTI 

This chapter is organized as follows. In Sec. I5.2I we overview the linearized Bekenstein 
model. In Sec. I5.3l we discuss the non-linear evolution of the fine structure constant using two 
different approximations. Our results are described in Sect. I5.4l and further discussion will be 
found in Sect. I5.5I 



Let us consider the linearized Bekenstein model (see for example Ref. H112I0 . In this model a 
neutral scalar field is non-minimally coupled to electromagnetism by means of the gauge kinetic 
function, Bp(ip)- The Lagrangian density is given by 



5.2 The Model 




(5.1) 



where 




(5.2) 



(53) 



and £ ot her is the Lagrangian density of the other fields. We assume that V((p) is a linear function 
of the field ip, namely, 

V{<p) = V(<p ) + ^- (<p - <p Q ) , (5.4) 
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where dV/dip is constant. The subscript '0' refers to the present time. Also, we consider that 
B F is a linear function of cp, 

B F {ip) = 1 - /% - ip ) • (5.5) 

Equivalence Principle tests impose constraints on the value of (3, namely \(3\ < 1CT 3 G 1 ^ 2 . We 
recall that B F (ip) acts as the effective dielectric permittivity which can be phenomenologically 
taken to be an arbitrary function of ip. 

Using the definition of the gauge kinetic function, Bp = a /a, we note that 

oi 

— = l + P(cp-<p ) (5.6) 
a 

up to linear order. Assuming that 

^<0, ^>0, (5.7) 
dip dip 

the results lead to a fine- structure constant that had a smaller value in the past than at the present 
time. The assumptions of linear functions of ip for V((p) and a ((f) are necessarily verified for 
some period of time around the present day. 

The equation of motion for the field cp is given by 



dV{<p) P 
dip 4 



Oip = + 9F^F^ (5.8) 



and we define 



-1^- — = - 7jF p m , (5.9) 

where ^ F represents the matter density fraction which contributes to the right hand side of Eq. 
(15.81) . For baryons 7^ is constrained to be very small since the electromagnetic corrections to the 
mass of protons and neutrons are of the order of 1 part in 10 4 [11 1311 . Also, since the variations 
of a are very small, 7^ is expected to be a nearly constant parameter. Bearing in mind that 

F^F^ = 2(B 2 - E 2 ) < 0, (5.10) 



we note that the last two terms in Eq. (|5.8I) have opposite signs. 
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In this chapter we will consider the linearized Bekenstein model with Eq. (I5.8I) written as 



□V5 = -/?7FPm, (5.H) 

where we have set dV/dtp = for simplicity. Of course a generalization to models with 
dV/dip 7^ would modify the background dynamics of a but would not otherwise change 
our conclusions. 

First, by using Eq. (15.111) we can compute the evolution of the background value of a with 
physical time, that is given approximately by 

a + 3Ha ~ — £,p m a , (5.12) 



where 

f = lF (3 2 < 6 x 10" 6 7 F . (5.13) 

Note that we have used the definition of the field B F in terms of a to write Eq. (|5.12l) . On 
the other hand, since we aim in this chapter to provide a description for the evolution of a in 
the non-linear regime by using a local static approximation, we linearize Eq. (15.111) to obtain a 
local static solution in a slightly perturbed Minkowski space 

V 2 a 

— ~£fym- (5-14) 

a 

Here, we neglected the <p term in Eq. (15.1 II) in order to obtain equation Eq. (15.141) . This is a 
valid approximation if | (p \ <C |V 2 <^|. 

Note that Eq. (15.121) implies that the time variation of a induced in one Hubble time is 

— ~ Z5p m aH- 2 , (5.15) 
a 

while Eq. (15.141) tell us that the spatial variation of a seeded by a perturbation of size L H^ 1 
is given by 

Z5 Pm aL 2 . (5.16) 

a 

Hence, we expect that the spatial variations of a generated on small cosmological scales (where 
non-linear gravitational clustering is currently taking place) will be small. 
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5.3 Non-Linear Evolution of a 

We will consider the evolution of a due to the growth of a uniform spherical matter inhomo- 
geneity with a final size r using two different approaches. 



5.3.1 Infinite Wavelength 

We will start by considering the evolution of two homogeneous and isotropic universes. This 
approach was used in It86*ll87l to study the evolution of a uniform (that is, infinite wavelength) 
perturbation in a, with the initial conditions set at some initial time ti deep in the matter era. 
Thus we name this the infinite wavelength approximation. 



5.3.1.1 The background universe 



The background is a spatially flat Friedman-Robertson- Walker universe containing matter, ra- 
diation and a cosmological constant. The Friedmann equation is 



H 2 = Hi 



0; 



+ Q r 



Ai 



(5.17) 



where VL mi , Q ri and are respectively the density parameter for matter, radiation and cos- 
mological constant at some initial time U. In order to ensure the flatness for the background 
universe, we assume that 

1 _ &"mi ~ &>ri — &>Ai = 0. (5.18) 

The acceleration of the universe is given by 



a 



— ^ + n ri (— 



-n 



Ai 



(5.19) 



Here for simplicity we include the contribution of the energy density of the field if in the value 
of thus neglecting the contribution of the kinetic contribution to energy density of the field 
Lf, which is nevertheless constrained to be small. Also, since any variation of the fine-structure 
constant from the epoch of nucleosynthesis onwards is expected to be very small ll45l l47ll we 
neglect its minor contribution to the evolution of the baryon density, included in VL m . 
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5.3.1.2 The perturbed universe 



We also consider a second universe, a spatially c 
with matter, radiation and a cosmological constani. 



osed Friedmann-Robertson- Walker universe 
3 with 



n 



hi 



1 Qrni ^ri ^Ai 



(5.20) 



This closed universe is obtained by perturbing the matter density of the background universe 
with a small infinite wavelenght perturbation. We take the initial time U to be deep into the 
radiation epoch in such way that 



H Pi ~ Hi, 

^Pri ^ri-j 

QpAi = ^Ai, 

Qpmi Q-mi AVL m i, 

^Pki AVl m i. 



(5.21) 
(5.22) 
(5.23) 
(5.24) 
(5.25) 



where the subscript P indicates that the quantity is related to the perturbed (closed) universe. 

In this approximation the values of a in the background and perturbed universes are com- 
puted using Eq. (15.121) so that 



a 



b,oo 



a 



b,oo 



3£l m H 2 



57T 



(5.26) 



where 



Aa a — OLi 



(5.27) 



a OLi 

and the subscripts b and oo indicate that a is calculated for the background and perturbed uni- 
verses respectively. 



'We have neglect the contribution of the scalar field <p for the background as well as for the closed universe but 
it can be included. This does not modify our conclusions. 
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5.3.2 Local Approximation 



We have introduced another approach in Ref. H which we have referred to as the local ap- 
proximation. In this approach, we consider a slightly perturbed Minkowski space and assume 
that the background value of the fine structure constant is fixed. Of course, this is a good ap- 
proximation only on length scales much smaller than the Hubble radius. Although we follow 
the evolution of a into the non-linear regime, where the density perturbations become much 
larger than unity, the fractional variations of a are always constrained to be small. Thus we 
name this the local approximation [|112] . 

We solve the Poisson equation (|5.14l) outside the spherical distribution of mass, M, with 
radius r, to determine the spatial variation of a as 



5a 
a 



a[r) — a 
a 



M 
Anr 



(5.28) 



identifying a with a(r = oo). Given that 



M 



Aix5p m r J 



(5.29) 



we have 



— = -l&Pmr 2 . 
a S 



(5.30) 



This calculation assumes that the spherical inhomogeneity has a radius significantly smaller 
than the horizon in order for Eq. (15.141) to be a valid approximation locally. It is straightforward 
to show that, in this limit, the condition 



d 2 (5a/a) 



dt 2 



5a 



a 



H 2 < 



d 2 (5a/a) 



dr 2 



(5.31) 



is required by self-consistency. 

Therefore Aa / a can be computed as we have 



Aa 

a 



Aa 



a 



+ 



5a 



a 



(5.32) 



If one considers the matter density computed for the background and perturbed universes in the 
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Figure 5.1: Evolution of the scale factor, a, as a function of physical time, t, for the background and 
perturbed universes (solid and dashed lines respectively) 



previous subsection, the local variation of a can be written as 



5a AH 2 n m ) P -(H 2 n m ) 2 
~ ~ e 8^ (aproF) 



(5.33) 



where the P represents the closed perturbed universe. We can already anticipate that these 
two different approaches for estimating spatial variations of a will produce very different re- 
sults since in the first approach there is no reference to the size of the fluctuation (an infinite 
approximation is considered). 



5.4 Numerical Results 

The initial values of the various cosmological parameters were chosen in such way that at 
present time we have = 0.29, VIa = 0.71 and f2 r0 = 8.4 x 10~ 5 for the background 
universe. Also, the value of Afi mi is such that the collapse of the perturbed universe occurs 
near the present epoch. This is clearly seen in Fig. 15.11 which plots the evolution of the scale 
factor a as as a function of the physical time, t, for the background and perturbed universes 
(solid and dashed lines respectively). Note that, according to the spherical collapse model, a 
uniform density perturbation will become infinite in a finite time. Clearly this is not realistic (at 
least for local fluctuations) since in practice the collapse is stopped by virialization. 
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Figure 5.2: Evolution of the fine structure constant, a, as a function of physical time, t, in the vicinity of 
spherical distribution of mass including the spatial variations of a calculated using the local approxi- 
mation with ro = 2Hq (solid line) and the infinite wavelength approximation (dashed line). The dotted 
line represents the background evolution of a. Note that £ < 6 x 10 -6 7i? is a very small number. 
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Figure 5.3: Evolution of (a — a&) /(<*& — (*{), in various models as a function of physical time, t. 
The dashed line represents the infinite wavelength approximation while the solid lines use the local 
approximation with ro = 2Hq , ro = Hq 1 and rg = Hq 1 /2 (top to bottom respectively). 
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We plot in Fig. I5.2I the b evolution of Aa/a obtained using the two different models. 
The background evolution is shown by the dotted line, while Aa/a for infinite wavelength 
and local approximation models are plotted with dashed and solid lines, respectively. The 
crucial difference among the two approximations is seen immediately by solving the equations 
of evolution: while the evolution of Aa/a in the local approximation depends on the final 
size of the perturbation, r , no reference to the perturbation size is made in the second case. 
Actually, for the latter model a uniform infinite wavelength is assumed. Hence, the latter model 
will produce a single result while the first method will produce results which depend explicitly 
on the final perturbation size, r . For example, the results shown in Fig. I5.2I for the local 
approximation are obtained by assuming r = IH^ 1 . If one works in the non-linear regime, 
this is clearly an unrealistic value of r . However, we have used it in order to get a sizable 
depart from the background Aa/a (note that in order to be self-consistent the model requires 
that r <C H^ 1 ) and the result can easily be rescaled for realistic values of r using (15.331) . A 
more realistic value of r is the radius of a typical cluster of galaxies (r ~ lOMpc). In this case 
the spatial variations of a calculated using the local approximation need to be rescaled down by 
roughly six orders of magnitude with respect to the unrealistic r ~ 2 Hq 1 case shown in Figs. 



For observing explicitly the behavior of the evolution of the fine-structure constant when 
smaller values of r are considered, that is, when it is assumed mass distributions into the non- 
linear regime we compute the results for (a — a&) / (c^ — 014) using the different models as a 
function of physical time, t. In Fig. 15.31 the dashed line represents the one result obtained by 
using the infinite wavelength approximation and the solid lines represent the results obtained 
by the local approximation with r = 2iJ ~ 1 , r = H^ 1 and r = H 1 /2 (top to bottom 
respectively). Although we have started by non-realistic values of the radius of the spherical 
distribution, it is clear that as we move towards smaller values of r the local approximation 
will give negligible variations of a, that is indeed the correct conclusion for realistic cases. 




5.5 Discussion 



It is interesting to discuss the evolution of 8a /a in the linear regime during the matter and 
radiation eras. In the linear regime the physical radius of the spherical overdensity grows as 
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r oc a, and consequently M varies as 

M oc Sp m a 3 = 5 m p m a 3 . (5.34) 

We conclude that the spherical mass distribution is approximately constant during the radiation 
era and grows proportionally to a during the matter era. 



Using the local approximation one has that according to the solution of the Poisson equation 
outside the overdensity, deep into the radiation era 

5a i 

— oc a , (5.35) 

a 



while deep into the matter era 



5a _ x 
— ~ constant . (5.36) 
a 



It is also important to point out that the amplitude of the a variation (for a given 5p m ) de- 
pends quadratically on the radius of the spherical overdensity (Aa/a oc r 2 ). Consequently, for 
r <C Hq 1 the spatial variations predicted by the local approximation will be very small. The 
background evolution of Aa/a in linearized Bekenstein models has been previously studied in 
HI 1211 and it is well known that Aa /a — > deep into the radiation era. 



This is confirmed by the results plotted in Fig. 15 .21 which clearly confirm the above discus- 
sion. In both approximations 

^(^0 (5.37) 

a V " J b 

deep into the radiation era. On the other hand 

Aa\ f Aa\ 

~ constant (5.38) 

b 



a ,, \ a 



deep into the matter era. Also, as expected the evolution of Aa/a in the infinite wavelength 
approximation only departs significantly from the background one near the present time (when 
the collapse of the perturbed universe occurs). 



146 Evolution of a in the Non-Linear Regime 



5.6 Conclusions 

We have presented and contrasted two different approximations for the non-linear evolution of 
a. While in the local approximation we assume that the size of the perturbation which gives 
rise to spatial variation of a is much smaller than the Hubble radius, H~ l , in the infinite wave- 
length approximation the opposite is required for self-consistency. Hence, it is not surprising 
that the results obtained using the two models are so different. The crucial question is which one 
provides the right answer when applied to cosmology. Clearly, the infinite wavelength approxi- 
mation is not a good approximation in general because we are not interested in an homogeneous 
perturbation with a size larger than the Hubble radius. Note that in this approximation not only 
the gradients are neglected but it is also assumed that the Hubble parameter, H, is the local one 
(not the global H). In addition, the fact that the matter density is homogeneously distributed 
inside a given spherical region does not imply that a will also be homogeneous within that 
region. 

Therefore we are convinced that local approximation provides the right answer since the 
non-linear effects are only expected to be important in this context on scales much smaller than 
Hq 1 (scales smaller than the typical galaxy cluster size). On such small scales, we therefore 
predict that the spatial variations of a generated in the simplest models should be too small to 
be of cosmological interest. 

These results confirm that it is difficult to generate significant large-scale spatial variations 
of a as it is shown in previous chapter for the special cases of cosmic strings and magnetic 
monopoles 03 13 even when we account for the evolution of the fine structure constant in the 
non-linear regime, when density perturbations become much larger than unity. This results were 
also confirmed by other authors in Ref. H109H - H111L 



Chapter 6 



Outlook 

In this thesis we have studied the cosmological consequences of topological defects, in partic- 
ular for the dark energy scenario and varying fundamental constants theories. In the following, 
we summarize our main results and conclusions and address future works. 

In Chapter O we have considered the possibility of a domain wall network as a source of 
dark energy to explain the recent acceleration of the universe. 

• We have presented an analytic model that describes very well the domain wall evolution. 
We have taken into account the curvature and the interactions of the walls with radiation 
and other particles and derived a strong bound on the curvature of the walls, which shows 
that viable candidate networks must be fine-tuned and non-standard. 

• We have studied several typical domain wall models in two spatial dimensions, and dis- 
cussed the conditions under which various types of defect junctions can exist.We have 
highlighted the situations where only stable Y- type junctions or X-type ones are formed. 
In addition we have shown the case where both types co-exist. 

• Using geometrical, topological and energy arguments we distinguished which features 
of a particular model are crucial to construct the best candidate potential for frustration. 
However, our results strongly lead us to a no frustration conjecture. In all our numerical 
simulations, including those of the ideal model, we have found a dynamical behavior 
consistent with a convergence towards a scaling solution and no evidence of any behavior 
that could eventually lead to frustration. 

• The results of 3D simulations clearly indicate that no frustrated network can be formed 
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dynamically out of random initial conditions that would mimic a realistic phase transi- 
tion. Therefore we have presented the most compelling evidence to date that domain wall 
networks can not be the dark energy. 

In this context, it is interesting to study the dependence of our results on the number of 
space-time dimensions. Another relevant problem is the role of the junctions in the dynamics 
of the network. We recall that the type of junctions formed will depend both on energetic 
considerations and on the topology of the minima of the potential in field space. In the future 
we intend to consider further examples of particular models and study this issue with analytic 
and numerical tools. 

In Chapter[3]we have investigated cosmic strings in the context of varying-a theories. 

• We have considered a generic gauge kinetic function and shown that there is a class of 
models of this type for which the classical Nielsen-Olesen vortex is still a valid solution. 
However, in general, the solutions of the local cosmic strings depart from the standard 
case. 

• We have shown that the Equivalence Principle constraints impose tight limits on the al- 
lowed variations of a on cosmological scales induced by cosmic string networks of this 
type. 

In Chapter |4] we have investigated the variations of the fine-structure constant in the vicinity of 
a static local magnetic monopole. 

• We have introduced the Bekenstein-type models in the Non-Abelian Yang-Mills theory 
and studied various models with varying-a and confirmed that despite the existence of a 
class of models for which the 't Hooft-Polyakov standard solution is still valid. 

• We have verified that, in general, the solutions of the local magnetic monopoles depart 
from the standard case. 

• We have shown that Equivalence Principle constraints impose tight limits on the vari- 
ations of a induced by magnetic monopoles which confirms the difficulty to generate 
significant large-scale spatial variation of the fine structure constant found in the previous 
section, even in the most favorable case where these variations are seeded by magnetic 
monopoles. 
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In Chapter \5\ we have investigated the evolution of the fine- structure constant in the non-linear 
regime. 

• We have presented and contrasted two different approximations for the non-linear evolu- 
tion of a. While in the local approximation we have assumed that the size of the perturba- 
tion which gives rise to spatial variation of a is much smaller than the Hubble radius, H* 1 , 
in the infinite wavelength approximation the opposite is required for self-consistency. 

• We have shown that the infinite wavelength approximation is not a good approximation 
in general because we are not interested in an homogeneous perturbation with a size 
larger than the Hubble radius. Note that in this approximation not only the gradients 
are neglected but it is also assumed that the Hubble parameter, H, is the local one (not 
the global H). In addition, the fact that the matter density is homogeneously distributed 
inside a given spherical region does not imply that a will also be homogeneous within 
that region. 

• We have concluded that local approximation provides the right answer since the non- 
linear effects are only expected to be important in this context on scales much smaller 
than Hq 1 (scales smaller than the typical galaxy cluster size). On such small scales, we 
therefore predict that the spatial variations of a generated in the simplest models should 
be too small to be of any cosmological interest. 

• We have presented results which confirm that it is difficult to generate significant large- 
scale spatial variations of a as it is shown in previous chapter for the special cases of 
cosmic strings and magnetic monopoles even when we account for the evolution of the 
fine structure constant in the non-linear regime, when density perturbations become much 
larger than unity. 

Overall, we have assumed that the spatial variations of a in the vicinity of cosmic strings, 
magnetic monopoles and compact objects are sourced only by the gauge kinetic function Bp 
coupled to the Maxwell term. A generalization of this model, we intend to investigate in the 
near future, includes different source of variation of the fine- structure constant due to a discrete 
spontaneous symmetry breaking which produces domain walls dividing the universe in patches 
with different values of a. Of course, the symmetry breaking scale must be sufficiently low 
in order for the domain walls to be cosmologically benign. This is an interesting model, in 
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which domain walls may be responsible for cosmologically relevant spatial variations of the 
fine- structure constant as well as other fundamental parameters. 



Appendix - Numerical Codes - Domain 
Wall Network Evolution 



The Main Code 

In this section we present the C code used for performing the field theory numerical simulation of domain 
wall networks evolution. This version carries out the evolution of ideal class of models for N = 20 and 
was optimized for the shared memory architeture of the COSMOS supercomputer were the simulations 

were performed. For a 512 3 box, this simulation took about 2 hours and 15 minutes on 128 processors. 

/* */ 

/* VT 26/10/06 vlO */ 

/* Compile with -DRAW for binary data output */ 

/* Compile with -DTEST for timings & fixed seed=0 */ 

/* */ 

#include <stdio.h> 
#include <math.h> 
#include <string.h> 
#include <stdlib.h> 
#include "walls.h" 
#ifdef _OPENMP_ 
#include <omp.h> 
#endif 

int main(int argc, char **argv) 
{ 

static double S0[diml][dim2][dim3][dim4][inicampos]; 
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static double S 1 [dim 1 ] [dim2] [dim3] [dim4] [inicampos] ; 
static double DS0[diml][dim2][dim3][dim4] [inicampos]; 
static double DSl[diml][dim2][dim3][dim4] [inicampos]; 
static double GS[diml][dim2][dim3][dim4] [inicampos]; 
double GT; 

static float GTf[diml][dim2][dim3]; 
double VS[inicampos+l]; 
double K[23] [23]; 
double R[inicampos+l],cons; 
int m,wq,qw; 

^integration variables inicialization*/ 

int i,ipl,iml; 

intjjpljml; 

int k,kpl,kml; 

int l,lpl,lml,NQ; 

int ncounter; 

/*evolution variables inicialization*/ 

double ct,dct,dx,dx2; 

double ctfinal; 

double at; 

double delta; 

double D1,D2; 

/*area and energy calculation inicialization*/ 

double V0t4; 

double miu; 

double A,Pi; 

double Lshi[inicampos]; 

double ek; 

double AREA, area; 

double difl,dif2,dif3,dif4; 

double VO, V; 

char doutfilme[MAXBUF], doutscaling[MAXBUF]; 
/* files*/ 

FILE *dados, *foutfilme, *scaling, *fcp, *fstep; 
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/*simulation parameters inicialization*/ 



int Nl= 


=diml; 


int N2= 


=dim2; 


int N3= 


=dim3; 


int N4= 


=dim4; 



double lamda; 
long idum; 
int skipcp = 0; 
char sstep[MAXBUF]; 
char scp[MAXBUF]; 
double nfieldsl; 
double epsilon=-0. 1 ; 
double omega0=5.0; 
int num,teste; 
double alfa=3.0; 
double beta0=0.0; 
double W0=10.0; 
double a,b,e; 

/*dynamical range inicialization*/ 
double zi; 
double tO,ti; 

int counterl, contar,D,c,s, cont, cont2; 

double expoentec,alpha; 

double eklim,ekaux,ekauxQ; 

int npontos,npontosQ; 

double assim=1.0; 

double r=0.5; 

int arg = (argc < 2 ? 1 : atoi(argv[l])); 
char nomeficheiro[20]; 
/* OpenMP */ 

double my_eklim, my_AREA, my_ek; 
int my_npontos; 
/* timing */ 

double tl,t2,t3,t4,t5, ta,tb,tc,td; 
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double csecond(); 
tl = csecond(); 
#ifdef TEST 
idum = 0; 
#else 

idum = time (NULL); 
#endif 

srand48(idum); 

/* get parameters from cmd line or default ini file */ 

if ( argv[l] != NULL && ((dados=fopen(argv[l],"r")) != NULL) ) { 

fscanf(dados,"%lf %lf %lf %lf %d %d %d %d n", 

& lim, & expoente,& ctinitial,& tampasso,& npassos,& nfields,& noutput,& nckpt); 

} else if ( (dados=fopen(inifile,"r")) != NULL ) 

{ 

fscanf(dados,"%lf %lf %lf %lf %d %d %d %d n", 

& lim, & expoente,& ctinitial,& tampasso,& npassos,& nfields,& noutput,& nckpt); 
} 

print_parameters() ; 
contar=0.0; 

dct=tampasso*ctinitial; 

ctfinal=ctinitial+npassos*dct; 

expoentec=expoente/( 1 -expoente) ; 

Pi=3.1415927; 

dx=ctinitial; 

dx2=dx*dx; 

ct=ctinitial; 

nfields 1 =nfields+ 1.0; 

ncounter=0; 

counter 1=0; 

lamda=0.3/(nfields+5.0); 

V0=0.145; 

c=0; 

^^^^^^^^^^^^^^^^^^ ^* ]\^xrnin3. ffitXtrix.^* h^^^^^^^^h^^^h^^^^^^^^^^^^^^^^^ 
K[l][l]=1.0; K[l][2]=0.0; K[l][3]=0.0; K[l][4]=0.0; K[l][5]=0.0; 
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K[l][6]=0.0; K[l][7]=0.0; K[l][8]=0.0; K[l][9]=0.0; K[l][10]=0.0; 
K[1][11]=0.0; K[l][12]=0.0; K[l][13]=0.0; K[l][14]=0.0; K[l][15]=0.0; 

K[l][16]=0.0; K[l][17]=0.0; K[l][18]=0.0; K[l][19]=0.0; K[l][20]=0.0; 

/***/ 

K[2][l]=- 1.0/2.0; K[2][2]=sqrt(3.0)/2.0; K[2][3]=0.0; K[2][4]=0.0; K[2][5]=0.0; 
K[2][6]=0.0; K[2][7]=0.0; K[2][8]=0.0; K[2][9]=0.0; K[2][10]=0.0; 
K[2][ll]=0.0; K[2][12]=0.0; K[2][13]=0.0; K[2][14]=0.0; K[2][15]=0.0; 

K[2][16]=0.0; K[2][17]=0.0; K[2][18]=0.0; K[2][19]=0.0; K[2][20]=0.0; 

/***/ 

K[3][l]=- 1.0/2.0; K[3][2]=-sqrt(3.0)/2.0; K[3][3]=0.0; K[3][4]=0.0; K[3][5]=0.0; 
K[3][6]=0.0; K[3][7]=0.0; K[3][8]=0.0; K[3][9]=0.0; K[3][10]=0.0; 
K[3][ll]=0.0; K[3][12]=0.0; K[3][13]=0.0; K[3][14]=0.0; K[3][15]=0.0; 

K[3][16]=0.0; K[3][17]=0.0; K[3][18]=0.0; K[3][19]=0.0; K[3][20]=0.0; 

/***/ 

K[4][l]=0.0; K[4][2]=0.0; K[4][3]=sqrt(2.0); K[4][4]=0.0; K[4][5]=0.0; 
K[4][6]=0.0; K[4][7]=0.0; K[4][8]=0.0; K[4][9]=0.0; K[4][10]=0.0; 
K[4][ll]=0.0; K[4][12]=0.0; K[4][13]=0.0; K[4][14]=0.0; K[4][15]=0.0; 

K[4][16]=0.0; K[4][17]=0.0; K[4][18]=0.0; K[4][19]=0.0; K[4][20]=0.0; 

/***/ 

K[5][l]=0.0; K[5][2]=0.0; K[5][3]=sqrt(2.0)/4.0; K[5][4]=sqrt(30.0)/4.0; K[5][5]=0.0; 
K[5][6]=0.0; K[5][7]=0.0; K[5][8]=0.0; K[5][9]=0.0; K[5][10]=0.0; 
K[5][ll]=0.0; K[5][12]=0.0; K[5][13]=0.0; K[5][14]=0.0; K[5][15]=0.0; 

K[5][16]=0.0; K[5][17]=0.0; K[5][18]=0.0; K[5][19]=0.0; K[5][20]=0.0; 

/***/ 

K[6][l]=0.0; K[6][2]=0.0; K[6][3]=sqrt(2.0)/4.0; K[6][4]=sqrt(30.0)/20.0; K[6][5]=3.0*sqrt(5.0)/5.0; 
K[6][6]=0.0; K[6][7]=0.0; K[6][8]=0.0; K[6][9]=0.0; K[6][10]=0.0; 
K[6][ll]=0.0; K[6][12]=0.0; K[6][13]=0.0; K[6][14]=0.0; K[6][15]=0.0; 

K[6][16]=0.0; K[6][17]=0.0; K[6][18]=0.0; K[6][19]=0.0; K[6][20]=0.0; 

/***/ 

K[7][l]=0.0; K[7][2]=0.0; K[7][3]=sqrt(2.0)/4.0; K[7][4]=sqrt(30.0)/20.0; K[7][5]=sqrt(5.0)/10.0; 
K[7][6]=sqrt(7.0)/2.0; K[7][7]=0.0; K[7][8]=0.0; K[7][9]=0.0; K[7][10]=0.0; 
K[7][ll]=0.0; K[7][12]=0.0; K[7][13]=0.0; K[7][14]=0.0; K[7][15]=0.0; 

K[7][16]=0.0; K[7][17]=0.0; K[7][18]=0.0; K[7][19]=0.0; K[7][20]=0.0; 

/***/ 
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K[8][l]=0.0; K[8][2]=0.0; K[8][3]=sqrt(2.0)/4.0; K[8][4]=sqrt(30.0)/20.0; K[8][5]=sqrt(5.0)/10.0; 
K[8][6]=sqrt(7.0)/14.0; K[8][7]=2.0*sqrt(21.0)/7.0; K[8][8]=0.0; K[8][9]=0.0; K[8][10]=0.0; 
K[8][ll]=0.0; K[8][12]=0.0; K[8][13]=0.0; K[8][14]=0.0; K[8][15]=0.0; 

K[8][16]=0.0; K[8][17]=0.0; K[8][18]=0.0; K[8][19]=0.0; K[8][20]=0.0; 

/***/ 

K[9][l]=0.0; K[9][2]=0.0; K[9][3]=sqrt(2.0)/4.0; K[9][4]=sqrt(30.0)/20.0; K[9][5]=sqrt(5.0)/10.0; 
K[9][6]=sqrt(7.0)/14.0; K[9][7]=sqrt(21.0)/28.0; K[9][8]=3.0*sqrt(3.0)/4.0; K[9][9]=0.0; K[9][10]=0.0; 
K[9][ll]=0.0; K[9][12]=0.0; K[9][13]=0.0; K[9][14]=0.0; K[9][15]=0.0; 

K[9][16]=0.0; K[9][17]=0.0; K[9][18]=0.0; K[9][19]=0.0; K[9][20]=0.0; 

/***/ 

K[10][l]=0.0; K[10][2]=0.0; K[10][3]=sqrt(2.0)/4.0; K[10][4]=sqrt(30.0)/20.0; K[10][5]=sqrt(5.0)/10.0; 
K[10][6]=sqrt(7.0)/14.0; K[10][7]=sqrt(21.0)/28.0; K[10][8]=sqrt(3.0)/12.0; 
K[10][9]=sqrt(15.0)/3.0; K[10][10]=0.0; 

K[10][ll]=0.0; K[10][12]=0.0; K[10][13]=0.0; K[10][14]=0.0; K[10][15]=0.0; 

K[10][16]=0.0; K[10][17]=0.0; K[10][18]=0.0; K[10][19]=0.0; K[10][20]=0.0; 

/***/ 

K[11][1]=0.0; K[ll][2]=0.0; K[ll][3]=sqrt(2.0)/4.0; K[ll][4]=sqrt(30.0)/20.0; K[ll][5]=sqrt(5.0)/10.0; 
K[ll][6]=sqrt(7.0)/14.0; K[ll][7]=sqrt(21.0)/28.0; K[ll][8]=sqrt(3.0)/12.0; K[l l][9]=sqrt(15.0)/30.0; 
K[ll][10]=sqrt(165.0)/10.0; 

K[11][11]=0.0; K[ll][12]=0.0; K[ll][13]=0.0; K[ll][14]=0.0; K[ll][15]=0.0; 

K[ll][16]=0.0; K[ll][17]=0.0; K[ll][18]=0.0; K[ll][19]=0.0; K[ll][20]=0.0; 

/***/ 

K[12][l]=0.0; K[12][2]=0.0; K[12][3]=sqrt(2.0)/4.0; K[12][4]=sqrt(30.0)/20.0; K[12][5]=sqrt(5.0)/10.0; 
K[12][6]=sqrt(7.0)/14.0; K[12][7]=sqrt(21.0)/28.0; K[12][8]=sqrt(3.0)/12.0; K[12][9]=sqrt(15.0)/30.0; 
K[12][10]=sqrt(165.0)/110.0; 

K[12][ll]=3.0*sqrt(22.0)/11.0;K[12][12]=0.0; K[12][13]=0.0; K[12][14]=0.0; K[12][15]=0.0; 

K[12][16]=0.0; K[12][17]=0.0; K[12][18]=0.0; K[12][19]=0.0; K[12][20]=0.0; 

/***/ 

K[13][l]=0.0; K[13][2]=0.0; K[13][3]=sqrt(2.0)/4.0; K[13][4]=sqrt(30.0)/20.0; K[13][5]=sqrt(5.0)/10.0; 
K[13][6]=sqrt(7.0)/14.0; K[13][7]=sqrt(21.0)/28.0; K[13][8]=sqrt(3.0)/12.0; K[13][9]=sqrt(15.0)/30.0; 
K[13][10]=sqrt(165.0)/110.0; 

K[13][ll]=sqrt(22.0)/44.0; K[13][12]=sqrt(26.0)/4.0; K[13][13]=0.0; K[13][14]=0.0; K[13][15]=0.0; 

K[13][16]=0.0; K[13][17]=0.0; K[13][18]=0.0; K[13][19]=0.0; K[13][20]=0.0; 

/***/ 
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K[14][l]=0.0; K[14][2]=0.0; K[14][3]=sqrt(2.0)/4.0; K[14][4]=sqrt(30.0)/20.0; K[14][5]=sqrt(5.0)/10.0; 
K[14][6]=sqrt(7.0)/14.0; K[14][7]=sqrt(21.0)/28.0; K[14][8]=sqrt(3.0)/12.0; K[14][9]=sqrt(15.0)/30.0; 
K[14][10]=sqrt(165.0)/110.0; 

K[14][ll]=sqrt(22.0)/44.0; K[14][12]=sqrt(26.0)/52.0; K[14][13]=sqrt(273.0)/13.0; K[14][14]=0.0; 
K[14][15]=0.0; 

K[14][16]=0.0; K[14][17]=0.0; K[14][18]=0.0; K[14][19]=0.0; K[14][20]=0.0; 

/***/ 

K[15][l]=0.0; K[15][2]=0.0; K[15][3]=sqrt(2.0)/4.0; K[15][4]=sqrt(30.0)/20.0; K[15][5]=sqrt(5.0)/10.0; 
K[15][6]=sqrt(7.0)/14.0; K[15][7]=sqrt(21.0)/28.0; K[15][8]=sqrt(3.0)/12.0; K[15][9]=sqrt(15.0)/30.0; 
K[15][10]=sqrt(165.0)/110.0; 

K[15][ll]=sqrt(22.0)/44.0; K[15][12]=sqrt(26.0)/52.0; K[15][13]=sqrt(273.0)/182.0; 
K[15][14]=3.0*sqrt(35.0)/14.0; K[15][15]=0.0; 

K[15][16]=0.0; K[15][17]=0.0; K[15][18]=0.0; K[15][19]=0.0; K[15][20]=0.0; 

/***/ 

K[16][l]=0.0; K[16][2]=0.0; K[16][3]=sqrt(2.0)/4.0; K[16][4]=sqrt(30.0)/20.0; K[16][5]=sqrt(5.0)/10.0; 
K[16][6]=sqrt(7.0)/14.0; K[16][7]=sqrt(21.0)/28.0; K[16][8]=sqrt(3.0)/12.0; K[16][9]=sqrt(15.0)/30.0; 
K[16][10]=sqrt(165.0)/110.0; 

K[16][ll]=sqrt(22.0)/44.0; K[16][12]=sqrt(26.0)/52.0; K[16][13]=sqrt(273.0)/182.0; 
K[16][14]=sqrt(35.0)/70.0; K[16][15]=2.0*sqrt(10.0)/5.0; 

K[16][16]=0.0; K[16][17]=0.0; K[16][18]=0.0; K[16][19]=0.0; K[16][20]=0.0; 

/***/ 

K[17][l]=0.0; K[17][2]=0.0; K[17][3]=sqrt(2.0)/4.0; K[17][4]=sqrt(30.0)/20.0; K[17][5]=sqrt(5.0)/10.0; 
K[17][6]=sqrt(7.0)/14.0; K[17][7]=sqrt(21.0)/28.0; K[17][8]=sqrt(3.0)/12.0; K[17][9]=sqrt(15.0)/30.0; 
K[17][10]=sqrt(165.0)/110.0; 

K[17][ll]=sqrt(22.0)/44.0; K[17][12]=sqrt(26.0)/52.0; K[17][13]=sqrt(273.0)/182.0; 
K[17][14]=sqrt(35.0)/70.0; K[17][15]=sqrt(10.0)/40.0; 

K[17][16]=sqrt(102.0)/8.0; K[17][17]=0.0; K[17][18]=0.0; K[17][19]=0.0; K[17][20]=0.0; 

/***/ 

K[18][l]=0.0; K[18][2]=0.0; K[18][3]=sqrt(2.0)/4.0; K[18][4]=sqrt(30.0)/20.0; K[18][5]=sqrt(5.0)/10.0; 
K[18][6]=sqrt(7.0)/14.0; K[18][7]=sqrt(21.0)/28.0; K[18][8]=sqrt(3.0)/12.0; K[18][9]=sqrt(15.0)/30.0; 
K[18][10]=sqrt(165.0)/110.0; 

K[18][ll]=sqrt(22.0)/44.0; K[18][12]=sqrt(26.0)/52.0; K[18][13]=sqrt(273.0)/182.0; 
K[18][14]=sqrt(35.0)/70.0; K[18][15]=sqrt(10.0)/40.0; 

K[18][16]=sqrt(102.0)/136.0; K[18][17]=3.0*sqrt(51.0)/17.0;K[17][18]=0.0; K[18][19]=0.0; 
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K[17][20]=0.0; 

/***/ 

K[19][l]=0.0; K[19][2]=0.0; K[19][3]=sqrt(2.0)/4.0; K[19][4]=sqrt(30.0)/20.0; K[19][5]=sqrt(5.0)/10.0; 
K[19][6]=sqrt(7.0)/14.0; K[19][7]=sqrt(21.0)/28.0; K[19][8]=sqrt(3.0)/12.0; K[19][9]=sqrt(15.0)/30.0; 
K[19][10]=sqrt(165.0)/110.0; 

K[19][ll]=sqrt(22.0)/44.0; K[19][12]=sqrt(26.0)/52.0; K[19][13]=sqrt(273.0)/182.0; 
K[19][14]=sqrt(35.0)/70.0; K[19][15]=sqrt(10.0)/40.0; 

K[19][16]=sqrt(102.0)/136.0; K[19][17]=sqrt(51.0)/102.0; K[19][18]=sqrt(57.0)/6.0; K[19][19]=0.0; 

K[19][20]=0.0; 

/***/ 

K[20][l]=0.0; K[20][2]=0.0; K[20][3]=sqrt(2.0)/4.0; K[20][4]=sqrt(30.0)/20.0; K[20][5]=sqrt(5.0)/10.0; 
K[20][6]=sqrt(7.0)/14.0; K[20][7]=sqrt(21.0)/28.0; K[20][8]=sqrt(3.0)/12.0; K[20][9]=sqrt(15.0)/30.0; 
K[20][10]=sqrt(165.0)/110.0; 

K[20][ll]=sqrt(22.0)/44.0; K[20][12]=sqrt(26.0)/52.0; K[20][13]=sqrt(273.0)/182.0; 
K[20][14]=sqrt(35.0)/70.0; K[20][15]=sqrt(10.0)/40.0; 

K[20][16]=sqrt(102.0)/136.0; K[20][17]=sqrt(51.0)/102.0; K[20][18]=sqrt(57.0)/114.0; 

K[20][19]=sqrt(570)/19.0; K[20][20]=0.0; 

/***/ 

K[21][l]=0.0; K[21][2]=0.0; K[21][3]=sqrt(2.0)/4.0; K[21][4]=sqrt(30.0)/20.0; K[21][5]=sqrt(5.0)/10.0; 
K[21][6]=sqrt(7.0)/14.0; K[21][7]=sqrt(21.0)/28.0; K[21][8]=sqrt(3.0)/12.0; K[21][9]=sqrt(15.0)/30.0; 
K[21][10]=sqrt(165.0)/110.0; 

K[21][ll]=sqrt(22.0)/44.0; K[21][12]=sqrt(26.0)/52.0; K[21][13]=sqrt(273.0)/182.0; 
K[21][14]=sqrt(35.0)/70.0; K[21][15]=sqrt(10.0)/40.0; 

K[21][16]=sqrt(102.0)/136.0; K[21][17]=sqrt(51.0)/102.0; K[21][18]=sqrt(57.0)/114.0; 
K[21][19]=sqrt(570)/380.0; K[21][20]=3.0*sqrt(70.0)/20.0; 
/* initialize in parallel to get first touch right - vt */ 
#pragma omp parallel for default (none) 

private(i,j,k,l,cont) shared(Nl,N2,N3,N4,nfields,S0,Sl,DS0,DSl,GS, GTf) 

for(i=0;i<Nl;i++){ 

for(j=0;j<N2;j++){ 

for(k=0;k<N3;k++){ 

GTf[i][j][k] = 0.0; 

for(l=0;l<N4;l++){ 

for (cont=0;cont<nfields;cont++) { 
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SO[i][j][k][l][cont]=0.0; 

DSO[i][j][k][l][cont]=0.0; 

Sl[i][j][k][l][cont]=0.0; 

DSl[i][j][k][l][cont]=0.0; 

GS[i][]][k][l][cont]=0.0; 

} } } } } 

^ miticil conditions 
sprintf( sstep, "%s-%d-%d", ckpt_prefix, MDIM, nfields); 
if ( (fstep = fopen(sstep, "r")) == NULL ) 
{ 

/* Init from random distribution */ 

for(i=0;i<Nl;i++){ 

for(j=0;j<N2;j++){ 

for(k=0;k<N3;k++){ 

for(l=0;l<N4;l++){ 

cons= drand48(); 

for (cont2=0;cont2<=(nfields) ;cont2++) { 

if ((cons>(cont2/(nfieldsl))) && (cons<((cont2+l)/nfieldsl))){ 

for (cont=0;cont<nfields;cont++){ 

DSO[i][j][k][l][cont]=0.0; 

S0[i] [j ] [k] [1] [cont]=K[cont2+ 1 ] [cont+ 1 ] ; 

}}}}}}} 
} else { 

/* printf("Have checkpoint n"); */ 
fscanf( fstep, "%d %lf", &ncounter, &ct ); 

sprintf( scp, "%s-%d-%d-%.2f.bin", ckpt_prefix, MDIM, nfields, ct); 

if ( (fcp = fopen( scp, "rb")) != NULL) 

{ 

fread(&S0[0][0][0][0][0], sizeof (double), nfields*Nl*N2*N3*N4, fcp); 
fread(&DS0[0][0][0][0][0], sizeof (double), nfields*Nl*N2*N3*N4, fcp); 
fclose(fcp); 
} else { 

printf("Error: Cannot open checkpoint data file %s n", scp); 
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exit(l); 
} 

/* adjust contar so output file numbering is correct */ 

if ( ncounter > ) 

contar=ncounter/noutput; 

fclose(fstep); 

printf("Restarting with ct = %g , ncounter = %d, contar = %d n", 
ct, ncounter, contar); 
skipcp = 1 ; 
} 

filename(doutscaling, doutscaling_prefix,nfields) ; 

at=pow(ctinitial,expoentec); 

/* timing */ 

t2 = csecond(); 

ta = t2 -tl; 

printf("Init: time= %5.2f sec n", ta); 

/* MAIN LOOP */ 

#pragma omp parallel default(none) 

private(i,j ,k,l,cont,m,ipl ,iml jpl jml ,kpl ,kml ,lp 1 ,1ml , 

GT,Lshi,R,VS,my_AREA,my_ek,my_eklim,V,ekaux) 

shared(Nl,N2,N3,N4,nfields,nfieldsl,GS,S0,dx2,delta,dct,lamda, 

Sl,DSl,K,lim,area,DS0,ek,eklim,ct,ctfinal,tl,t2,t3,t4, 

ta,tb,tc,td,ncounter,alfa,expoentec,contar,doutfilme, 

foutfilme,at,scaling,omegaO,VO, noutput, nckpt, fstep, 

sstep, scp, fcp, skipcp, doutscaling, GTf, doutfilme_prefix, 

ckpt_prefix) 

{ 

while (ct < ctfinal){ 
#pragma omp master 
{ 

#ifdef TEST 
ta=t2-tl; 
tb=t4-t2; 
tc=0; 
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td=t4-tl; 

printf ("Step:%3d ta=%5.2f tb=%5.2f tot=%5.2f ", 
ncounter, ta,tb,td); 
tl = csecond(); 
#endif 

delta=0.5*alfa*dct/ct*expoentec; ek=0.0; 

area=0.0; 

eklim=0.0; 

} 

/* end master */ 
#pragma omp barrier 
my_ek = 0; 
my_eklim = 0.0; 
my_AREA=0; 
#pragma omp for 
for(i=0;i<Nl;i++){ 
for(j=0;j<N2;j++){ 
for(k=0;k<N3;k++){ 
for(l=0;l<N4;l++){ 
ipl=(i+l) %N1; 
iml=(i-l+Nl) % Nl; 
jpl=(j+l)%N2; 
jml=(j-l+N2) % N2; 
kpl=(k+l) % N3; 
kml=(k-l+N3) % N3; 
lpl=(l+l) % N4; 
lml=(l-l+N4) % N4; 
V=0.0; 
GT=0.0; 

/* GS from SO; GT is sum of GS (over cont) */ 
for (cont=0;cont<nfields;cont++){ 

GS[i]|j][k][l][cont]=0.5*(pow((0.5*(S0[ipl]|j][k][l][cont]-S0[iml]|j][k][l][cont])),2.0) 
+pow((0.5*(S0[i] [jpl] M [1] [cont]-S0[i] [jml] [k] [1] [cont])),2.0) 
+pow((0.5*(S0[i] [j] [kp 1] [1] [cont]-S0[i] [j] [kml] [1] [cont])),2.0) 
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+pow((0.5*(S0[i][j][k][lpl][cont]-S0[i][j][k][lml][cont])),2.0)); 

GT+=GS[i][j][k][l][cont]; 

} 

GTf[i][j][k] = (float) GT; 

/* Lshi from SO (over cont) */ 

for (cont=0;cont<nfields;cont++){ 

Lshi[cont]=SO[ipl] [j] [k] [1] [cont]+S0[iml] [j] [k] [1] [cont] 

-2.0*S0[i] [j] [k] [1] [cont]+S0[i] [jp 1] [k] [1] [cont] 

+S0[i] [jml] [k] [1] [cont]-2.0*S0[i] [j] [k] [1] [cont] 

+S0[i] [j] [kpl] [1] [cont]+S0[i] [j] [kml] [1] [cont] 

-2.0*S0[i] [j] [k] [1] [cont]+S0[i] [j] [k] [lpl ] [cont] 

+S0[i] [j] [k] [lm 1 ] [cont] -2.0*S0[i] [j] [k] [1] [cont] ; 

Lshi[cont]/=dx2; 

} 

/* R from SO (over m over cont) */ 
for(m=0;m<nfields 1 ;m++) { 
R[m]=0; 

for (cont=0;cont<nflelds;cont++){ 
R[m]+=pow(SO[i][j][k][l][cont] - K[m+l][cont+l],2.0); 
VS[cont]=0; 
} } 

/* VS from R (over m over cont) */ 

for(m=0;m<nfields 1 ;m++) { 

for (cont=0;cont<nfields;cont++){ 

VS[cont]+=lamda*3.0*(R[m]-3.0)*(R[m]-1.0)*(S0[i][j][k][l][cont]-K[m+l][cont+l]); 
} } 

ekaux = 0.0; 

/* DS1 from DSO & VS & Lshi, SI from SO & DS1 (over cont)*/ 

for (cont=0;cont<nfields;cont++){ 

DS 1 [i] [j] [k] [1] [cont]=((l .0-delta)*DS0[i] [j] [k] [1] [cont] 

+dct* (Lshi [cont] - VS [cont] ))/( 1 .0+delta) ; 

SI [i] [j] [k] [1] [cont]=S0[i] [j] [k] [1] [cont]+dct*DS 1 [i] [j] [k] [1] [cont] ; 
ekaux+=DS 1 [i] [j] [k] [1] [cont] *DS 1 [i] [j] [k] [1] [cont] ; 



} 

for(m=0;m<nfields 1 ;m++) { 

V+=lamda*0.5*R[m]*(R[m]-3.0)*(R[m]-3.0); 

} 

if (V >= lim){ 
my_eklim += ekaux/V; 
my_AREA+=1.0; 
} } } } } 

#pragma omp for 

for(i=0;i<Nl;i++){ 

for(j=0;j<N2;j++){ 

for(k=0;k<N3;k++){ 

for(l=0;l<N4;l++){ 

for (cont=0;cont<nfields;cont++){ 

S0[i][j][k] [1] [cont] =S 1 [i] [j] [k] [1] [cont] ; 

DS0[i] [j] [k] [1] [cont]=DS 1 [i] [j] [k] [1] [cont] ; 

my_ek+=DS 1 [i] [j] [k] [1] [cont] *DS 1 [i] [j] [k] [1] [com 

} } } } } 

#pragma omp atomic 
area += my_AREA; 
#pragma omp atomic 
ek += my_ek; 
#pragma omp atomic 
eklim += my_eklim; 
#pragma omp barrier 
#pragma omp master 
{ 

#ifdef TEST 
t2 = csecond(); 
#endif 
ct+=dct; 

at=pow(ct,expoentec) ; 
ncounter+=l; 

if ( (scaling=fopen(doutscaling,"a")) != NULL ) 
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{ 

fprintf(scaling,"%e %e %e %e n'\ct,area/(diml*dim2*dim3*dim4),eklim/area,ek/((area)*omega0*V0)); 

fclose(scaling); 

} 

if (((ncounter-1) % noutput )==0.0) { 

printf("(out)"); 

contar+=l; 

filename(doutfilme,doutfilme_prefix,nfields*100+contar); 
if ( (foutfilme=fopen(doutfilme,"w")) != NULL ) { 
#ifdef RAW 

fwrite(GTf, sizeof(float), N1*N2*N3, foutfflme); 
#else 

for(i=0;i<Nl;i++){ 
for(j=0;j<N2;j++){ 
for(k=0;k<N3;k++){ 
fprintf(foutfilme,"%f n",GTf[i] [j] [k]); 

} } } 
#endif 

fclose(foutfilme); 
} } 

/* checkpoint every nckpt - except if just restored */ 

if ( skipcp == 1 ) 

skipcp = 0; 

else 

{ 

if ( ncounter>0 && nckpt>0 && (ncounter % nckpt == ) ) { 

/* printf("Checkpointing at ct=%g, ncounter=%d ", ct, ncounter);*/ 

printf("(ckpt)"); 

sprintf( scp, "%s-%d-%d-%.2f.bin", ckpt_prefix, MDIM, nfields, ct); 

if ( (fcp = fopen( scp, "w")) != NULL) 

{ 

fwrite(&S0[0][0][0][0][0], sizeof (double), nfields*Nl*N2*N3*N4, fcp); 
fwrite(&DS0[0][0][0][0][0], sizeof (double), nfields*Nl*N2*N3*N4, fcp); 
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fclose(fcp); 
} else { 

printf("Error: Cannot open file %s n", scp); 

exit(l); 

} 

sprintf( sstep, "%s-%d-%d", ckpt_prefix, MDIM, nfields); 

if ( (fstep = fopen(sstep, "w")) != NULL ) 

{ 

fprintf( fstep, "%d %lf", ncounter, ct ); 

fclose(fstep); 

} else { 

printf("Error: Cannot open file %s n", sstep); 

exit(l); 

} } } 

printf(" n"); 
#ifdef TEST 
t4 = csecond(); 
#endif 
} 

#pragma omp barrier 
} 

/* end loop over ct */ 
} 

*/ end parallel*/ 
} 

Simulation Parameters & default values 

include "walls.h" 
/* Filenames */ 

char inifile[MAXBUF] = "dados.txt" ; 
char doutfilme_prefix[MAXBUF] = "box3D" ; 
char doutscaling_prefix[MAXBUF]= "scaling3D" ; 
char ckpt_prefix[MAXBUF] = "ckpt3D" ; 
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/* Default values */ 

double lim = 0.0290 ; 

double expoente = 0.6666666667; 

double ctinitial = 1.0; 

double tampasso = 0.25 ; 

int npassos = 1020 ; 

/* number of steps */ 

int nfields = inicampos ; 

/* needs to be <= inicampos*/ 

int nckpt = 255 ; 

/* frequency to checkpoint */ 

int noutput = 1 ; 

/* frequency to output field array */ 

void print_parameters(void ) 

{ 

printf("# ####################################### n n "); 

printf("Simulation parameters: n n"); 

printf(" dim = %d n ", MDIM); 

printf(" lim = %lf n", lim); 

printf(" exp = %lf n", expoente); 

printf(" ct_ini = %lf n", ctinitial); 

printf(" tampasso = %lf n", tampasso); 

printf(" steps = %d n ", npassos); 

printf(" fields = %d n", nfields); 

printf(" n output = %d n", noutput); 

printf(" n ckpt = %d n", nckpt); 

printf(" n"); 

printf("# ####################################### n n "); 
} 



Filename Input 

void filename(outputname, inputname, outnum) 
char outputname[]; 
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char inputname[]; 
int outnum; 
{ 

int dig4, dig3, dig2, digl; 

int inull; 

register int jout; 

dig4 = outnum / 1000; 

dig3 = (outnum - 1000 * dig4) / 100; 

dig2 = (outnum - 1000 * dig4 - 100 * dig3) / 10; 

digl = outnum - 1000 * dig4 - 100 * dig3 - 10 * dig2; 

inull = 1 ; 

while (((inputname[inull - 1] !='')& & (inputname [inull - 1] != ' 0') & & (inull < 26))) 
inull = inull + 1 ; 

for (jout = 1; jout <= inull - 1; jout++) 
outputname[jout - 1] = inputname [jout - 1]; 
outputname[inull - 1] = '.'; 
outputname[inull + 1 - 1] = 48 + dig4; 
outputname[inull + 2 - 1] = 48 + dig3; 
outputname[inull + 3 - 1] = 48 + dig2; 
outputname[inull + 4 - 1] = 48 + digl; 
outputname[inull + 5 - 1] =' 0'; 
} 

Dimensions and Number of Fields 

# define MDIM512 

# define diml MDIM 

# define dim2 MDIM 

# define dim3 MDIM 

# define dim4 1 
/* preprocessor */ 

# define MAXFUB 100 
/* function definitions */ 

extern void filename ( char [], char [], int); 
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extern double csecond ( void); 
extern void print_parameter( void); 
/* filename */ 

extern char inirile[MAXFUB]; 

extern char doutfilme_prefix[MAXFUB]; 

extern char doutscaling_prefix[MAXFUB]; 

extern char ckpt_prefix[MAXFUB]; 

/*default values*/ 

extern double lim; 

extern double expoente; 

extern double ctinitial; 

extern double tampasso; 

extern double npassos; 

extern double nfields; 

extern double nckpt; 

extern double noutput; 
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